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Abstract 

The  purpose  of  this  research  is  to  determine  if  acoustical  effects  on  a  coaxial  shear 
injection  jet  flow  can  be  modeled  utilizing  molecular  dynamic  simulation.  Molecular 
dynamic  simulations  model  flows  as  a  group  of  interacting  particles.  The  flow  in  this 
research  was  simulated  using  nitrogen  molecules.  The  initial  task  involved  achieving 
effective  geometry  for  a  simulated  coaxial  jet.  The  coaxial  jet  geometry  was  driven  by 
the  desire  for  simulations  to  operate  in  the  continuum  regime,  which  requires  very  low 
Knudsen  numbers.  Three  outer  to  inner  jet  ratios  of  0.0,  1.0,  and  6.0  were  examined  with 
the  inner  jet  velocity  maintained  constant  at  50  m/s. 

Velocity  profiles  in  the  injector  need  to  be  controlled  in  order  to  validate  the  continuum 
flow.  Acoustic  interference  is  introduced  into  the  simulation,  and  mixing  and  density 
profiles  provide  valuable  infonnation  into  the  how  the  flow  is  affected  by  the  acoustic 
interference.  Radial  density  profiles  also  provide  infonnation  about  the  shape  of  the  jet 
with  and  without  acoustic  interference  as  it  exits  the  injector.  The  effects  of  acoustic 
interference  for  most  cases  showed  good  agreement  with  previous  experimental  data. 
Results  showed  good  validation  of  the  simulation  and  warrants  more  in-depth  study. 
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MODELING  ACOUSTIC  EFFECTS  ON  SHEAR  COAXIAL  JET  FLOW 
UTILIZING  MOLECULAR  DYNAMIC  SIMULATION 


I.  Background 


Chapter  Overview 

This  chapter  presents  a  brief  look  into  the  motivation  for  this  research.  A  description  of 
the  problem  of  acoustic  interference  on  coaxial  jet  flow  is  discussed  along  with  some 
effects  of  the  acoustic  interference.  The  research  objectives  are  also  reviewed.  These 
objectives  will  provide  guidelines  for  a  successful  research  effort.  The  results  presented 
later  in  the  document  also  reflect  the  objectives,  and  provide  an  answer  to  some  of  the 
problems  presented  in  our  problem  description. 

Problem  Description 

Jet  flow  has  been  a  topic  of  research  for  decades.  Coaxial  jet  flows  adds  a  higher  level  of 
complexity  to  jet  study.  Adding  the  study  of  acoustic  effects  brings  an  even  greater  level 
of  complexity  to  jet  flow  research.  Historically,  research  involving  coaxial  jets  and 
acoustic  effects  has  taken  place  in  a  physical  laboratory  with  expensive  equipment  and 
computational  tools  used  to  capture  the  desired  relationships  between  a  coaxial  jet’s  core 
and  outer  flows.  Similar  tools  can  be  used  to  measure  the  inherent  acoustics  present  in 
operating  engines,  and  the  effects  of  these  acoustics  on  the  flow  issuing  out  from  the  jet. 
This  experimentation  can  become  costly  especially  when  the  need  for  multiple  test 
conditions  exists.  One  solution  points  to  the  possibility  of  modeling  these  same  jet  flows 
using  a  computer  simulation. 


Previous  research  has  seen  molecular  dynamics  used  to  simulate  supercritical  fluid  flow 
(Branam,  2005).  The  first  challenge  becomes  adding  a  coaxial  component  to  the 
simulation.  Once  the  coaxial  flow  is  modeled  accurately,  an  operating  geometry  must  be 
determined  to  minimize  the  computationally  intensive  workload.  The  ideal  ratio  of  outer 
to  inner  jet  area  needs  to  be  determined  for  the  molecular  dynamic  simulation. 

Once  the  molecular  dynamic  simulation  is  capable  of  both  coaxial  flow  and  acoustic 
interference,  certain  criteria  must  be  met  to  detennine  if  the  molecular  dynamic 
simulation  may  be  used  to  accurately  model  the  problem  of  acoustic  interference  in  a 
coaxial  injection  jet. 

Objectives  of  Research 

Utilizing  molecular  dynamics  to  model  the  acoustic  effects  on  a  shear-coaxial  jet  flow 
requires  some  initial  work  and  experimentation.  The  velocity  profile  for  our  inner  and 
outer  tubes  needs  to  be  controllable.  The  minimum  ratio  of  tube  length  to  tube  diameter 
necessary  to  achieve  the  predetennined  theoretical  velocity  profile  for  the  simulations 
must  be  found.  There  is  a  range  of  tube  length  to  diameter  ratios  providing  the  desired 
velocity  profile  but  due  to  computational  and  memory  restrictions,  the  minimum  tube 
length  to  diameter  ratio  providing  the  required  velocity  curve  is  preferable. 

Adding  shear-coaxial  jet  flow  and  acoustic  interference  capabilities  to  the  molecular 
dynamics  tool  allows  this  research  to  achieve  its  objectives.  Shear-coaxial  jet  flow 
requires  new  geometry  considerations  and  increases  the  needed  computing  power.  The 
acoustic  interference  models  are  to  be  verified  by  examining  the  effects  of  jet  mixing 
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with  the  interference  turned  on  and  off.  Adding  acoustic  interference  should  increase 
mixing  once  the  jet  injects  into  the  quiescent  chamber. 

Once  the  coaxial  and  acoustic  capabilities  are  available,  the  objective  becomes  comparing 
flow  fields  with  and  without  acoustic  interference.  Density  profiles,  which  give  insight 
to  the  shape  of  a  jet,  are  examined.  Jets  without  acoustic  interference  will  show  steeper 
density  curves  along  with  higher  densities  in  the  center  of  the  jet.  Once  acoustic 
interference  is  added,  the  jet  is  expected  to  have  a  flatter  density  profile  while  becoming 
more  diffuse  and  dissipating  quicker. 

The  simulation  explores  a  range  of  outer  to  inner  jet  velocities.  The  velocity  profiles  for 
the  inner  and  outer  jets  will  provide  data  that  will  aid  in  validating  the  flow  issuing  into 
the  chamber.  Figure  1  gives  a  graphical  representation  of  the  problem  description. 


I  Acoustic  Interference  1 


Figure  1.  Problem  Description 
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Thesis  Scope 

This  document  contains  four  main  chapters.  Chapter  two  presents  the  literature  review 
including  the  theory  behind  the  multiple  disciplines  taken  into  account  for  the  research. 
The  chapter  looks  into  coaxial  jet  theory,  including  velocity  and  diameter  differentials 
between  the  inner  and  outer  jets.  Acoustic  interference  is  discussed  as  it  relates  to  jet 
streams.  A  discussion  of  molecular  level  acoustic  scaling  and  nanoflows  is  also 
necessary  and  is  included  in  this  chapter. 

Chapter  three  covers  the  research  methodology,  discussing  molecular  dynamic  simulation 
in  detail,  considering  this  technique  is  the  backbone  of  the  research.  The  chapter  reviews 
force  potential  models,  wall  modeling,  and  numerical  integration  techniques.  Techniques 
aiding  in  decreasing  computation  times  such  as  parallel  computing  are  also  discussed. 

The  method  used  to  size  the  coaxial  jet  presents  itself  here  as  well.  The  chapter  also 
looks  at  the  various  simulated  cases  studied  for  this  research,  such  as  test  cases  flowing 
with  and  without  a  chamber.  The  research  shows  comparisons  for  coaxial  jets  with 
varying  velocity  differences  between  the  inner  and  outer  jets.  The  research  highlights  the 
effects  of  acoustics  on  these  cases. 

Chapter  four  presents  the  results  of  the  research  while  chapter  five  brings  conclusions  and 
recommendations  for  bettering  and  continuing  the  research  effort. 
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II.  Literature  Review 


Chapter  Overview 

This  chapter  includes  a  review  of  the  many  disciplines  incorporated  into  this  research 
effort.  The  literature  review  begins  with  a  simple  discussion  of  narrow  axisymmetric  jets 
followed  by  a  discussion  of  laminar  and  turbulent  jets  and  transition  theory.  Section 
three  delves  into  coaxial  jet  theory.  A  coaxial  flow  field  is  broken  down  to  increase 
understanding  of  the  effects  of  a  coaxial  jet  on  the  flow  cone  especially  with  respect  to 
the  mixing  of  the  flow  and  regions  of  change  in  the  flow.  The  effects  of  coaxial  jet 
geometry  on  a  flow  is  reviewed,  along  with  previous  research  relating  the  effects  of 
velocity  differences  between  the  inner  and  outer  jet  flow.  A  review  of  work  with  varied 
coaxial  jet  velocity  ratios  is  also  highlighted.  Section  four  reviews  acoustic  theory  and 
focuses  mainly  on  plane  waves.  This  section  will  also  develop  general  building  blocks  of 
acoustics  and  sound  wave  theory  equations  needed  to  develop  acoustic  interference  on  a 
molecular  scale.  Plane  waves  create  the  acoustic  interference  in  the  molecular 
simulations.  A  brief  discussion  of  nano  flow  theory  appears  at  the  close  of  the  chapter. 

Narrow  Axisymmetric  Jets 

If  a  round,  narrow  jet  issues  from  a  small  injector  with  sufficient  momentum,  it  will 
remain  narrow  and  grow  slowly  (White,  1974:350).  The  jet  changes  longitudinally  as 
dldx  and  radially  as  d/  dr ,  the  larger  of  the  changes.  Figure  2  represents  the  round, 
narrow  jet  issuing  from  an  injector. 
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Figure  2.  Round  Jet,  Free  Shear  Flow  (Wilcox,  2006:60) 

Squire  discovered  an  exact  solution  for  this  case.  Squire’s  stream  function  solution  for 
small  9  reduces  to 


lim  f{0)  =  lim- 

<?->0  v  ’  0->O  1 


2  sin2  6 


£ 


+  a-  cos  9  1  +  ^/4 


where  £  = 


261 


(1) 


Equation  2  shows  the  boundary  layer  solution  for  the  narrow  jet  with  the  origin  existing 
at  the  jet  exit. 

0  =  tan  (y  /  x)  <§ cl  (2) 

For  the  reduced  stream  function,  9  =  y/x  and  a  is  a  constant  less  than  0.05  and  is  directly 
related  to  the  dimensionless  jet  momentum.  The  longitudinal  velocity  follows  from  the 
definition  of  the  stream  function  (equation  3)  where  J  is  the  jet  momentum. 


■9)= 


3  J 


8  K/UX 


i+L 


(3) 


The  term  inside  of  the  parentheses  in  the  longitudinal  velocity  equation  represents  the 
shape  of  the  jet  while  the  jet  momentum  affects  the  jet  width.  This  term  also  represents 
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the  shape  of  a  turbulent  round  jet.  Equation  4  shows  mass  rate  of  flow  across  any  section 
of  the  flow  ( m  )  increases  linearly. 


oo 

m  =  /?  J  u  2  k  ycly  =  8  nfix 
0 


(4) 


Laminar  and  Turbulent  Jets 

Laminar  flow  has  a  very  low  resistance  to  high  Re  (White,  1974:385).  The  critical  Re  is 
defined  as  the  Re  where  a  laminar  flow  becomes  turbulent.  This  number  exists  on  the 
order  of  1,000  for  tube  flow,  not  a  large  value.  Due  to  this  low  critical  Re,  a  laminar  flow 
is  the  exception  to  the  rule  and  is  rarely  seen.  Due  to  its  low  resistance  to  Re,  a  laminar 
flow  is  in  essence  an  unstable  flow,  which  with  small  disturbances  moves  into  the 
turbulent  regime.  The  same  results  are  clearly  apparent  in  physical  experimentation 
(White,  1974:385-387).  If  laminar  flow  is  the  exception,  then  turbulence  must  exist  for 
the  majority  of  flows.  A  shear  laminar  jet  with  no  wall  or  pressure  gradients  can  be 
considered  relatively  unstable.  Most  boundary-layer  flows  are  stable  for  very  small 
disturbances  but  a  larger  disturbance  drives  the  flow  towards  instability.  This  transition 
from  laminar  to  turbulent  is  explained  using  small-disturbance  stability  theory 
(White,  1974:387).  Small-disturbance  theory  analysis  begins  with  a  set  of  governing 
equations  with  a  known  solution.  A  small  disturbance  is  added  to  this  solution.  After 
some  algebra,  the  small  disturbance  terms  encapsulate  the  original  solutions  and  become 
the  only  terms  left.  The  solution  is  simplified  by  assuming  that  the  flow  is  in  one 
direction.  The  simplification  reduces  the  equations  to  a  set  of  ordinary  homogeneous 
differential  equations  for  the  small  variation. 
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Following  the  small-disturbance  stability  theory  by  considering  the  Navier-Stokes 
equations  for  two  variables,  the  two  initial  equations  become 

V-F  =  0  (5) 

—  =  --Wp  +  W2V  (6) 

Dt  p 

A  solution  for  laminar  flow  is  set  as  VQ=(U,  V,  W)  =  Va  (x,  t )  and  pa  (x,  t )  and  the  small 
disturbance  is  defined  as  v(x,t)  =  (u,v,w)  and  p(x,t ) .  “All  disturbances  have  a  wave 
number  a,  a  propagation  speed  c,  and  a  frequency  co=ac  (White,  1974:392).”  Assuming 
that  the  flow  is  in  one-direction  helps  decrease  the  number  of  ordinary  differential 
equations  to  one.  The  result  is  the  Orr-Summerfeld  relation  found  in  equation  7,  a  fourth 
order  linear  homogeneous  equation  which  holds  the  key  to  infinitesimal  laminar-flow 
stability  (White:  1974:396). 

(U  -  c)(v"  -  a2v)  -  U"v  + — (v""  -  2a  V  +  a  V)  =  0  (7) 

The  boundary  conditions  for  the  Orr-Summerfeld  relationship  have  the  small 
disturbances  vanishing  at  infinity  and  at  any  walls  representing  a  no  slip  condition.  The 
boundary  conditions  that  apply  to  our  simulations  are  those  for  the  free  shear  flow,  which 
are 

v(±oo)  =  v'(±oo)  =  0  (8) 

Due  to  the  order  and  type  of  relationship,  the  Orr-Sommefeld  relationship  can  be 
approached  as  an  eigen  function  and  eigen  value  problem.  Orr-Summerfeld  simply  aids 
in  determining  the  stability  of  any  parallel  flow  along  a  single  direction  which  can  be 
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found  from  a  two-dimensional  analysis  for  the  effective  basic  flow  U0(y)  in  that  direction 
(White,  1974:396). 

Coaxial  Jets 

“There  is  a  frequent  need  to  realize  a  unifonn  mixture  from  two  initially  segregated 
streams  in  many  practical  instances  (Villermaux,  1998:807).”  A  coaxial  jet  is  one  of  the 
solutions  to  the  problem  of  mixing  two  streams.  Coaxial  jet  flow  can  occur  in  one  of  two 
ways.  First,  the  injector  may  itself  be  three-dimensional  as  in  a  rectangular  or  elliptical 
shape.  There  could  also  be  multiple  injectors  exhausting  parallel  to  the  main  flow  field 
operating  closely  enough  so  their  mixing  zones  interfere  (Schetz,  1980: 137).  The  second 
definition  of  a  coaxial  jet  is  the  one  that  is  of  importance  for  this  research  effort  and  will 
be  the  focus  of  the  research  work.  The  coaxial  jet  studied  has  an  injector  surrounded  by  a 
solid  tube  with  another  circular  injector  surrounding  the  solid  tube.  The  coaxial  jet  will 
inject  into  a  quiescent  chamber. 

The  coaxial  jet  flow  field  is  characterized  by  three  zones,  the  potential  core  (PC)  region, 
the  characteristic  decay  (CD)  region  and  the  axisymmetric  decay  (AD)  region.  The  PC 
region  is  where  “the  mixing  initiated  at  the  boundaries  has  not  yet  permeated  through  the 
entire  flow  field  leaving  a  region  of  unifonn  velocity  close  to  the  jet  exit  velocity  (Schetz, 
1980:137).”  The  potential  core  region  contains  two  potential  cores  and  two  distinct 
mixing  regions.  The  two  potential  cores  correspond  to  the  inner  and  outer  cores  issuing 
into  the  flow  stream.  The  first  mixing  region  occurs  between  the  boundary  of  the  inner 
stream  and  the  outer  stream.  The  second  mixing  region  occurs  between  the  outer  stream 
and  the  ambient  atmosphere  into  which  the  stream  issues.  The  CD  region,  which  is  the 
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region  of  axis  velocity  decay,  is  dependent  on  the  configuration  of  the  coaxial  jet.  The 
next  is  the  AD  region  ‘where  the  entire  flow  is  found  to  approach  symmetry  (Schetz, 
1980:137).’  This  region  is  where  the  flow  has  fully  mixed  and  has  reached  a 
homogeneous,  self-sustaining  flow  stream  that  becomes  identical  to  the  symmetrical  free 
jet  (Champagne  and  Wygnanski,  1971 : 1449).  As  the  stream  moves  into  the  AD  region, 
the  velocity  of  the  jet  becomes  easier  to  discern  as  it  now  depends  only  on  its  position 
and  total  jet  momentum  at  a  point.  Other  research  (Au  and  Ko,  1987:427)  has  referred  to 
these  three  zones  in  simpler  terms  as  the  initial  merging  zone,  the  intermediate  merging 
zone  and  the  fully  merged  zone.  Figure  3  shows  a  representation  of  the  coaxial  jet  flow 
flow  field  and  the  three  zones  within  the  flow  field. 


Potential  Core,  x„ 


Characteristic  Decay - 


.  Axisymmetric  Decay. 


Figure  3.  Characterization  of  Coaxial  Jet  Flow  Field  (Branam,2005) 
Champagne  and  Wygnanski  (1971 : 1445-1464)  discussed  the  idea  that  turbulent  flow  at 
some  point  becomes  fully  developed  and  self-preserving  and  the  velocity  at  a  point  is 
only  dependent  on  the  position  at  that  point  and  the  total  jet  momentum.  At  this  point, 
the  distribution  of  the  mean  velocity  and  the  turbulent  intensity  can  be  expressed  as  a 
single  variable  77  =  ytlh  where  y  is  the  radial  component  and  lh  is  the  length  scale.  The 
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scales  of  the  velocity  and  the  length  are  defined  as  U refDo(x  -  xh)  1  and  (x  -  xh ) ,  where 
xh  is  the  distance  of  the  virtual  origin  from  the  nozzle  exit  plane.  Reference  velocity  Uref 


can  be  calculated  using  total  jet  momentum  flux.  For  a  flow  ratio  U0/UiUp  to  10,  Uref  can 
be  defined  as  follows,  where  Jr  is  the  total  jet  momentum: 


(9) 


The  mean  velocity  and  turbulent  intensity  are  respectively  as  follows: 

—  =  :0 
Uref  x-xo  U  J 


(10) 


For  cases  where  U0  =  U,  or  U0  =  0,  the  mean  velocity  and  the  turbulent  intensity  can  be 
simplified.  For  U0/Ui  ranging  from  0.0  to  1.0,  the  inner  core  length  will  be  somewhat 
longer,  possibly  resulting  in  a  lower  shear  between  inner  and  outer  jets.  As  U0/Ui 
increases,  the  inner  potential  core  decreases.  Larger  U0/Ui  values  also  lead  to  higher 
turbulent  intensities  in  the  initial  regions  of  the  flow  (Champagne  and  Wygnanski, 
1971:1452). 

As  the  flow  issues  from  the  coaxial  jet,  changes  to  the  mean  velocity  ratio  do  not  affect 
the  outer  potential  core  length.  The  outer  potential  core  length  is  in  essence  a  constant, 
feeling  little  to  no  effect  from  changes  to  the  inner  jet  velocity  (Ko  and  Au,  1985:21 1- 
232;  Au  and  Ko,  1987:427-443).  Experimental  research  has  shown  the  width  of  each 
core  decreases  as  the  flow  moves  downstream  (Champagne  and  Wygnanski,  1971:1445- 
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1464).  The  ratio  of  annular  to  core  area  affects  the  pressure  differential  bending  the  jet 
inward  or  outward,  in  effect  shortening  or  lengthening  the  inner  core. 

Acoustic  Standing  Waves 

When  discussing  acoustics  in  a  traditional  sense,  the  focus  moves  towards  treating  fluid 
flow  as  a  continuum  using  macroscopic  quantities  such  as  density,  velocity  and  pressure 
as  dependent  variables  (Danforth-Hanford,  O’Connor,  Long,  and  Anderson,  2005).  A 
major  decision  must  be  made  between  two  methods  when  choosing  to  simulate  acoustics 
in  a  flow  field.  The  deciding  factor  has  been  shown  to  be  the  Knudsen  number  (Kn), 
where  Kn  is  defined  as  the  mean  free  path  divided  by  the  wavelength. 

Kn  =  -  (12) 

L 

There  are  separate  considerations  for  low  versus  high  Kn.  Sound  propagation  at  low  Kn 
relates  to  low  frequencies  and  atmospheric  conditions  while  sound  propagation  for  high 
Kn  relates  to  higher  frequencies  and  more  dilute  gases  (Danforth-Hanford,  et.  al,  2005). 
Danforth  and  Long  (2004:1948-1955)  have  shown  absorption  and  speed  of  sound  relies 
heavily  on  Kn.  The  continuum  model  of  fluid  flow  utilizes  the  Navier-Stokes  or  Euler 
equations  and  is  a  more  traditional  method  for  simulating  fluid  but  it  is  not  very  useful  for 
Kn  >  0. 1 .  The  particle  method,  while  more  complex,  provides  a  wider  operating 
spectrum  for  Kn  and  offers  the  ability  to  work  with  a  wider  range  of  fluid  flows.  Particle 
methods  are  necessary  for  Kn  >0.1,  which  can  occur  during  shockwaves,  high 
frequencies,  rarified  gas  or  other  flows  (Danforth  and  Long,  2004:1948). 

Looking  at  acoustics  in  the  audio  spectrum,  the  smallest  wavelength  possible  is 
approximately  2  cm  (Ludwig,  undated,  A).  A  significantly  smaller  wavelength  is  needed 
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in  order  to  bring  the  intended  level  of  interference  to  our  molecular  simulation.  A  cube 
of  air  with  a  width  approximately  one  thousand  times  smaller  than  the  2  cm  wavelength 
can  contain  over  2x10* 1  molecules  (Ludwig,  undated,  A).  Following  these  assumptions,  a 
cube  with  a  width  of  approximately  200A,  which  is  in  the  realm  of  interest  for  molecular 
dynamic  simulations,  could  contain  over  2xl07  molecules. 

The  molecular  dynamic  tool  models  the  acoustic  interference  in  the  simulations  in  two 
ways.  One  of  the  methods  will  have  the  interference  modeled  as  a  pressure  exerted  by 
the  gas  molecules  due  to  the  collisions  of  these  molecules  against  a  surface.  The 
“pressure  is  equal  to  the  change  in  momentum  caused  by  collisions,  per  unit  time  and  per 
unit  area  (Ludwig,  undated,  A).”  Total  momentum  flux  results  from  the  derivation  of  the 
momentum  flux  in  a  y-direction  across  the  face  of  the  cube  parallel  to  the  x-z  plane. 
Depending  on  the  interaction  of  the  acoustic  pressure  with  the  desired  surface,  the  flow 
across  the  remaining  four  walls  of  the  cube  may  need  to  be  considered  (most  of  the 
interference  in  the  research  simulations  takes  place  against  one  plane  of  the  chamber). 
Equations  13  and  14  can  be  combined  and  divided  through  by  AV  to  achieve  equation  15 
reflecting  the  law  of  conservation  of  momentum  (Ludwig,  undated,  B). 


a  r 


-  —  axavaz  \jig-m-s 


^(puxu  )  +  -^(puvuz) 
ox  oz 


axavaz  \jtg-m-s  ’’J 


dt 


(pu)  =  -V(pa2)-puV -u-(u •  V)(pu)  kg-m-s  2j 


The  pressure  term  in  the  preceding  equations  are  the  pa2  terms  where  P  =  pa2. 


(13) 

(14) 

(15) 
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Equation  16  expresses  the  Euler  equation,  an  expansion  from  the  momentum  equations. 
The  differential  equation  for  pressure  originates  from  the  energy  flux  equations  (equation 
17). 

~^~  =  — -V(pcr2)-(i7-V)i7  [m-s-2]  (16) 


d(p&2) 

dt 


7  2^7  - 

—  pa  V  ■  u 
5 


(17) 


The  pressure  differential  equations  and  the  Euler  equations  can  be  broken  into  two  main 
terms.  The  first  is  the  drift  velocities  reflecting  the  first-order  changes  directly  caused  by 
the  average  molecular  velocity  u  .  The  second  part  can  be  referred  to  as  the  diffusion 
tenns,  which  reflect  the  effects  of  molecular  diffusion  (Ludwig,  undated,  C).  For  almost 
all  cases,  the  diffusion  tenn  is  much  larger  than  the  drift  velocity  term.  Therefore,  the 
drift  velocity  term  drops  out  of  the  Euler  and  pressure  equations  leaving  the  following 
equations  (Ludwig,  undated,  C): 


du 

dt 


dp 

dt 


7 

5 


p,p:y  ■  w 


(18) 

(19) 


These  approximated  equations  represent  the  acoustic  equations. 

The  basis  of  acoustic  interference  still  falls  out  from  the  simpler  wave  fonnulas.  A 
review  of  the  following  equations  aids  in  building  towards  acoustic  pressure  amplitude 
determination.  The  following  are  basic  sound  wave  equations  where  f  is  frequency,  T  is 
the  oscillatory  period,  and  L  is  the  wavelength. 
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2  7t 

co  =  2nf  =  —  (angular  frequency) 


2  n 

k  =  —  (wave  number) 


ry  _  2nLf 
k  2  jz 


=  Lf  =  v  (speed  of  sound  through  medium) 


(20) 


(21) 


(22) 


Calculating  wave  number  and  angular  frequency  leads  into  the  pure  sine  wave  equation 

y  =  Asm[kx-cot^  (23) 

After  attaining  the  wave  attributes,  the  pressure  amplitude  of  the  wave  is  calculated  using 
the  following  equation  where  p0  is  the  density  of  the  medium  the  wave  travels  through. 

P„  =  P0covsso  (24) 


The  displacement  amplitude  (A)  is  directly  determined  from  the  specified  decibel  level 
by  using  equation  25  from  kinetic  theory  and  the  definition  of  the  speed  of  sound.  In  this 
expression,  co  is  the  angular  frequency  and  a  is  the  speed  of  sound  through  the  modeled 
medium. 


A  =  —  (25) 

pcoa 

Quantifying  the  pressure  amplitude  of  sound  was  necessary  to  create  an  acoustic  sound 
wave  for  simulation.  The  pressure  amplitude  determines  the  size  of  the  sound  signal. 
Pressure  amplitude  is  the  measure  of  the  size  of  the  variation  of  atmospheric  pressure 
caused  by  the  sound,  or  to  what  level  the  sound  pressure  disturbs  the  surrounding 
quiescent  air  (“Decibel  Scales”,  undated).  The  units  for  pressure  amplitude  for  sound  are 
decibels.  The  equation  for  sound  pressure  level  (SPL)  is 
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p 


(26) 


Lp  =  20  log 

\“o 

LP  represents  the  sound  pressure  level  in  dB,  P  represents  the  pressure  amplitude  of  the 
sound  disturbance,  and  Po  represents  a  reference  pressure  (2xl0‘  N/m  ).  It  is  important 
to  note  that  for  some  of  the  loudest  sound  disturbances,  air  pressure  varies  less  than  1%. 


Nanoflows 

Earlier,  Knudsen  number  was  defined  as  mean  free  path  divided  by  wavelength,  but  Kn 
can  also  be  equated  to 


Kn  = 


yn  M 
V  2  Re 


(27) 


which  represents  a  length  over  which  very  large  variations  of  a  macroscopic  quantity  can 
take  place  (Karniadakis,  Beskok,  and  Aluru,  2005:5 1-166).  For  rarified  gases,  Kn  plays  a 
large  role  in  the  flow  treatment.  Kn  in  the  ranges  of  0.01  to  0. 1  exhibit  slip  flow 
characteristics.  Transition  flow  characteristics  are  exhibited  when  Kn  spans  the  range  of 
0.1  to  10.  The  fluid  can  be  considered  a  continuum  flow  for  Kn  <0.01  (Karniadakis, 
Beskok,  and  Aluru,  2005:5 1-166).  The  molecular  dynamic  simulations  will  operate  in 
the  continuum  flow  regime. 

Utilizing  a  tangential  momentum  flux  relationship,  the  slip  velocity  is 


us =^h+(1-cr.) 


UA  +CTvUw~\ 


Velocity  profile  boundary  conditions  for  micro  flows  can  follow  two  models.  Model  A  is 
a  second  order  equation  while  Model  B  is  a  simpler  first  order  equation  with  a  general 
slip  coefficient,  b,  included  in  the  equation.  Models  A  and  B  are  both  non- 
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dimensionalized  with  respect  to  length  and  velocity  scale  and  are  shown  in  equations  29 
and  30  respectively  (Karniadakis,  et.  al,  2005:51-166). 


Us-Uw  = 


2-cr 


cr 


Kn 


f  du^ 

v  on  j 


+  - 


Kn 2 


dru 

V  dn  j 


(29) 


Us  -  Uw  = 


2  -cr„ 


Kn 


1  -  bKn 


f  dU} 
\  dn  j 


(30) 


For  flows  in  a  channel,  model  B  will  work  well  in  the  center  of  the  channel  but  will  give 
slight  deviations  near  the  walls.  These  observed  deviations  are  at  a  maximum  at  the  slip 
location.  Model  A  will  give  larger  errors  towards  the  center  of  the  channel.  The  velocity 
distribution  can  also  be  non-dimensionalized  as 


U *  (y ,  Kn)  =  U  (x,  y )  /U(x)  = 


\h) 


Kn 
h  1  -  bKn 


+  ^  + 


|  Kn 
6  1  -  bKn 


(31) 


where  b  =  0  represents  the  Maxwell  equation  of  the  first  order  and  b  =  - 1  represents  the 
second  order  boundary  condition  within  the  slip  regime  only.  The  velocity  boundary 
condition  models  and  non-dimensionalized  distribution  will  affect  the  Knudsen  layer,  a 
sub  layer  between  the  viscous  boundary  layer  and  the  wall  within  approximately  one 
mean  free  path. 
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III.  Methodology 


Chapter  Overview 

This  chapter  will  focus  on  simulation  set-up.  The  first  section  will  take  a  detailed  look  at 
the  molecular  dynamic  theory  and  inputs  used  for  each  simulation.  Molecular  dynamics 
topics  include  but  are  not  limited  to  soft  spheres,  Leonard-Jones  force  potential,  periodic 
boundary  conditions  and  molecular  wall  conditions.  The  following  section  will  provide  a 
breakdown  of  the  coaxial  jet  sizing  and  jet  flow  for  the  simulations.  The  last  section  will 
detail  the  acoustic  methods  used  in  the  simulations  to  provide  the  necessary  interference 
allowing  for  a  comparison  with  the  physical  experimentation  provided  in  previous 
research  (Davis,  2006). 

Molecular  Dynamics 

Molecular  dynamics  techniques  numerically  solve  the  N-body  problem  of  classical 
mechanics  (Haile,  1992:15).  Simply  put,  molecular  dynamics  provides  the  techniques  to 
solve  the  second-order  equation  of  motion 

m.r.  =  /  (32) 

l  l  J  l 

where  lnfis  the  atom’s  mass  (Allen  and  Tildesley,  1987:71).  Molecular  dynamic 
simulation  offers  the  ability  to  look  at  single  molecules  and  their  interaction  with  larger 
scale  phenomena.  Molecular  dynamic  simulation  can  be  used  to  understand  the  effects  of 
molecule  excitation  and  its  effects  on  the  flow  around  it  or  how  individual  molecules 
form  in  a  turbulent  flow.  The  simulation  process  on  a  molecular  level  can  be  broken  into 
three  steps.  The  first  step  requires  the  modeling  of  a  single  molecule  specific  to  the  type 
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of  simulation  and  flow  composition  required.  The  second  step  requires  the  modeling  of  a 
flow  consisting  of  many  of  the  modeled  molecules  from  step  one.  The  final  step  moves 
to  the  analysis  of  the  flow  data. 

Molecular  dynamic  particles  are  modeled  in  one  of  two  ways.  A  simulation  can  contain 
atoms  modeled  as  hard-spheres.  In  this  instance,  as  molecules  move  around  in  the  flow 
they  begin  to  collide,  reacting  similar  to  the  balls  in  a  game  of  billiards.  As  the  molecules 
hit,  there  is  no  deformation  of  the  molecules  and  the  collision  provides  an  immediate 
repulsion  of  the  molecules.  Soft-sphere  atoms  move  through  the  flow  like  hard  spheres 
but  the  soft  sphere  atom  interactions  evolve  differently.  Soft-spheres  are  modeled  using 
the  Lennard-Jones  potential  and  take  into  account  the  longer  range,  attractive  dispersion 
forces  along  with  the  short  range  repulsive  overlap  forces  maintaining  the  integrity  of  the 
molecules  themselves  (Haile,  1992:188-190).  This  simulation  models  molecules  as  soft- 
spheres  using  the  Lennard-Jones  potential. 

Molecular  dynamics  depends  on  the  interactions  of  N  atoms  in  a  flow,  where  N  could  be 
well  over  a  million  (1E6)  atoms.  The  possibility  of  more  than  (1E6)"  potential  atom 
interactions  would  prove  to  be  too  computationally  intensive  and  could  cause  simulations 
to  run  for  months  without  completion.  Instead,  the  interaction  energy  of  N  atoms  is  taken 
as  “the  sum  of  isolated  two-body  contributions”  (Haile,  1992:188-210).  Equations  33 
and  34  respectively  express  the  Lennard-Jones  potential  (u(r))  and  force  (F(r))  models 
using  two  body  interactions  for  soft-spheres,  where  a  represents  the  distance  to  the  zero 
in  u(r),  and  s  represents  the  energy  at  the  minimum  in  u(r). 

/  V2  f  \6 

<7  [Qi 

u{r)  =  4s  -  -  -  (33) 

\r )  yr )  _ 
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(34) 


Limiting  the  interaction  distances  among  the  atoms  decreases  computational  time  further. 
If  an  atom  sits  past  a  certain  distance  from  another,  the  two  atoms  have  no  significant 
interaction.  So  in  essence,  a  minimum  distance  is  selected  for  interaction  of  atoms 
causing  the  Lennard-Jones  potential  equation  to  be  expanded  to  equation  35  where  rc  is 
the  minimum  interaction  distance  for  the  atoms: 


The  simulation  flow  uses  nitrogen  as  the  modeled  molecule.  The  wall  fluid  potential 
interaction  (s/k)  for  nitrogen  is  36.5  K  and  the  location  of  the  well  between  two  particles 
(a)  is  3.293  A.  These  parameters  are  used  in  the  Lennard-Jones  force  routine  (Branam, 
2005:25). 


Verlet  Algorithm  and  RATTLE  Constraints 

The  Verlet  algorithm  solves  the  equations  of  motion  for  molecular  dynamic  simulation. 
The  method  is  a  direct  solution  of  the  second-order  equations  of  motion  shown  in 
equation  32.  The  Verlet  algorithm  method  maintains  very  good  energy  conservation 
properties  for  fairly  long  time  steps  and  is  easily  programmable.  The  Verlet  algorithm 
only  needs  to  account  for  the  position  r  (t),  accelerations  a  (t),  and  the  position  of  the 
previous  time  step  r  (t-5t)  (Allen  and  Tildesley,  1982:78-82).  The  velocities  themselves 
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are  not  present  due  to  cancellations  using  Taylor  expansions  and  are  not  necessary  to 
solve  for  trajectories. 

The  Verlet  algorithm  requires  two  previous  positions  to  estimate  the  present  position. 
The  algorithm  is  not  a  self-starting  system  and  requires  more  than  just  a  starting  position 
and  velocity  in  order  to  proceed.  Some  backward  method  such  as  a  backward  Euler 
method  must  be  used  in  order  to  get  r  (-5t)  (Haile,  1992: 158-159).  The  method  begins 
with  equations  36  and  37,  which  are  added  to  produce  the  final  Verlet  algorithm  for 
positions  which  is  seen  in  equation  38. 


fit* 

r  (7  + Jt)  =  r  (t)  +  St  •  v  (t)  +  •  a(f)  + ... 

(36) 

r[t-8t)  =  r  +  — 

(37) 

+  Jt)  =  2r  (t)~r  (t-  Jt)  +  St2  + 

(38) 

Even  though  it  contains  no  third  order  term,  the  equation  is  of  the  third  order.  The  final 
tenn  in  equation  38  expresses  the  local  truncation  error  which  is  of  the  order  (5t)4  (Haile, 
1992:158).  Even  though  the  velocity  is  not  needed  for  the  Verlet  algorithm  it  can  found 
by  using  equation  39. 


r  {t  +  dt^-r  {t  - dt~) 

2 Jt 


(39) 


Choosing  the  correct  step  size  is  also  important  when  using  the  Verlet  algorithm.  Using 
the  simulation  of  liquid  argon  near  the  triple  point  as  an  example,  RMS  fluctuations  were 
found  to  be  on  the  order  of  0.01  percent  of  potential  well  depth  for  a  time  step  5t«  10' 14 
seconds  and  the  fluctuations  increased  to  0.2  percent  for  5t«4xl0’14  seconds  (Allen  and 
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Tildesley,  1987:79-80).  The  simulation  time  steps  were  kept  at  10'14  seconds  or  10 
femtoseconds  and  checked  for  validity  against  analytical  solutions. 

The  basic  Verlet  algorithm  handles  velocities  in  a  very  awkward  manner.  Numerical  error 
can  be  introduced  when  the  basic  Verlet  algorithm  is  utilized.  The  velocity  Verlet 
algorithm  provides  solutions  to  the  shortcomings  of  the  basic  Verlet  algorithm.  The 
velocity  Verlet  algorithm  stores  positions,  velocities,  and  accelerations  all  at  the  same 
time  t,  which  minimizes  any  errors  due  to  round-off  (Allen  and  Tildesley,  1987:81)  The 
velocity  Verlet  algorithm  is  still  a  two  step  integration  and  only  introduces  truncation 
errors  of  the  fourth  order  (Branam,  2005:25-26).  First,  equation  40  is  used  to  calculate 
position  at  time  t+5t.  A  mid-step  velocity  is  calculated  using  equation  41.  At  this  point, 
forces  and  accelerations  are  calculated  at  the  next  time  step  and  the  new  velocity  at  step 
t+5t  is  calculated  from  equation  42.  The  numerical  stability,  convenience,  and  simplicity 
of  the  velocity  Verlet  algorithm  make  it  the  best  algorithm  option  presently  available 
(Allen  and  Tildesley,  1987:82).  Branam  (2005:25-27)  shows  the  velocity  Verlet 
algorithm  coupled  with  the  Lennard-Jones  potential  provides  more  accuracy  with  respect 
to  energy  conservation  when  compared  to  other  higher  order  integration  techniques. 
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Molecular  dynamic  simulations  require  constraints  be  applied  to  the  system  bonds  or  they 
can  be  modeled  as  discrete  potentials  separately.  While  many  of  an  atom’s  degrees  of 
freedom  change  under  the  molecular  forces  acting  on  the  system,  the  simulation  must 
maintain  the  integrity  of  an  atom’s  bond.  Algorithms  exist  to  offer  computationally 
efficient  bond  constraints.  An  approach  to  solving  constraints  involves  a  set  of  unknown 
multipliers  representing  the  magnitude  of  forces  along  the  bonds.  These  forces  are 
required  to  keep  the  bond  lengths  constant.  The  equations  of  motion  are  then  solved  for 
one  time  step  with  no  constraint  forces  present.  The  magnitude  of  the  constraint  forces 
can  then  be  determined  and  used  to  correct  the  atomic  positions  (Allen  and  Tildesley, 
1987:92-98).  The  RATTLE  algorithm  handles  the  equations  of  motion  when  internal 
constraints  are  present  in  the  system.  RATTLE  can  be  considered  a  variation  of  the 
Verlet  algorithm  with  internal  constraints  in  which  both  the  position  and  velocity  at  time 
step  t  are  used  to  calculate  the  position  and  velocity  at  time  step  t+5t  (Andersen,  1983). 
RATTLE  algorithm  allows  the  simulation  to  constrain  other  lengths  such  as  internal  bond 
angles  when  the  flow  problem  is  not  concerned  with  internal  molecular  modes  (Branam, 
2005:27). 

Periodic  Boundary  Conditions  and  Neighbor  List 

Periodic  boundary  conditions  model  a  repeating  pattern  in  molecular  dynamics  to 
simulate  infinite  dimensions  in  a  given  direction  when  modeling  small  systems  of  atoms. 
A  cell  created  that  is  smaller  than  the  simulation  volume  and  holds  N  number  of  atoms 
sets  periodic  boundary  conditions.  This  cell  becomes  the  primary  cell  and  is  considered  a 
small  portion  of  the  full  volume.  Exact  replica  cells  surround  the  primary  cell.  These 
cells  are  considered  image  cells,  as  they  are  in  fact  mirror  images  of  the  primary  cell. 
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These  image  cells  surround  the  primary  cell  and  hold  the  same  N  number  of  atoms.  The 
atoms  in  the  image  cells  move  as  copies  of  the  atoms  in  the  primary  cell.  Atoms  are  free 
to  pass  from  cell  to  cell  through  the  boundaries  but  each  cell  maintains  the  same  number 
of  molecules  because  as  a  molecule  leaves  one  cell,  a  similar  atom  is  entering  the  cell 
from  the  opposite  face.  The  computations  are  simplified,  because  the  only  molecules  we 
need  to  track  are  the  ones  in  the  primary  cell.  Using  the  example  of  a  simple  cubic  cell, 
all  other  molecule  positions  are  a  result  from  the  position  of  an  atom  in  the  primary  cell, 
length  of  the  cell  and  cell  transition  vector  a  as  seen  in  the  following  figure  (Haile, 
1992:81-86). 


Figure  4.  Application  of  periodic  boundary  conditions  in  two  dimensions  for  simple 
cubic  cells.  As  an  atom  leaves  the  primary  cell,  one  can  be  seen  entering  the  primary  cell 
from  the  opposite  face.  The  mirror  cells  mimic  the  atom  movements  taking  place  in  the 

primary  cell.  (Haile,  1992) 
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Other  molecules  located  within  a  certain  range  of  the  molecule  affect  the  interactions  of  a 
single  molecule  inside  a  simulation.  The  molecules  within  this  range  are  held  in  a 
variable  called  a  neighbor  list.  The  neighbor  list  is  another  mechanism  used  in  molecular 
dynamics  to  decrease  the  computational  time.  A  neighbor  list  limits  the  amount  of 
reactions  needing  consideration  by  one  molecule.  Each  molecule  located  inside  the 
neighbor  list  affects  the  other  molecules  inside  that  list.  The  reactive  or  attractive  force 
of  any  molecule  outside  of  the  neighbor  list  is  ignored  due  to  the  low  level  of  force 
potential.  The  list  is  not  a  stagnant  input  and  requires  updating  in  relation  to  the  position 
change  5r  of  the  molecules.  It  is  an  inverse  relationship,  as  5r  increases  the  neighbor  list 
refresh  time  decreases.  Typically,  neighbor  lists  will  be  updated  every  10-20  time  steps 
for  a  relatively  small  5r  (Rapaport,  1995). 

Molecular  Dynamics  Simulation  Development 

Branam  (2005)  created  and  first  used  the  molecular  dynamics  computer  simulation  code 
used  in  this  research.  The  code  modeled  supercritical  fluids  in  a  molecular  dynamic 
setting.  The  original  code  examined  the  operation  of  a  single  jet.  For  this  research,  the 
code  includes  the  option  for  studying  a  coaxial  jet  flow.  The  coaxial  jet  coding  allows  for 
varied  geometries  based  on  the  needs  of  the  user  and  the  extent  of  the  computer  resources 
available.  Coaxial  velocities  may  also  be  varied  as  desired  by  the  user.  An  acoustic 
interference  element  has  been  added  to  the  code  in  order  to  study  the  effects  on  the  jet 
flow.  Both  single  jet  and  coaxial  jet  flow  simulations  can  employ  the  acoustic 
interference  effects.  The  acoustic  interference  approach  appears  later  in  this  chapter. 
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Molecular  Dynamic  Simulation  Set-up 

The  molecular  dynamic  simulation  set-up  consisted  of  concentric  circular  jets  issuing  into 
a  rectangular  box-shaped  chamber  (Figure  1).  The  initial  research  consideration  was 
deciding  the  physical  dimensions  of  the  simulations  along  with  inner  and  outer  flow 
velocities.  The  theoretical  flow  velocity  used  for  this  effort  was  50  m/s  (Branam,  2005). 
The  simulations  were  run  using  Poiseuille  flow  to  force  particles  from  the  tube  into  the 
chamber  at  the  desired  speed.  A  theoretical  velocity  curve  was  determined  for  the 
pressure  driven  flow  at  the  velocity  of  50  m/s.  The  first  step  in  the  experimentation  is  to 
determine  the  minimum  length  over  diameter,  L/D  for  the  tube  in  order  to  achieve  the 
theoretical  velocity  curve  previously  detennined.  The  minimum  L/D  will  aid  in  keeping 
the  processing  time  down  while  providing  the  necessary  environment  for  valid  molecular 
dynamic  simulations.  For  this  initial  experimentation,  the  tube  radius  was  set  at  20 A  to 
keep  the  computation  time  to  a  minimum.  The  simulation  was  run  neglecting  the  chamber 
for  this  case.  The  coxial  tube  was  modeled  as  an  infinite  tube  in  order  to  match  the 
velocity  profde.  The  three  L/D  values  considered  at  this  tube  radius  were  1,  5,  and  10. 

The  force  multiplier  gives  the  simulation  the  indication  of  how  hard  to  push  the  flow. 
Using  this  information,  two  simulations  were  run  for  each  of  the  three  L/D  values.  Each 
L/D  value  was  simulated  with  a  negligible  force  multiplier  and  with  the  full  force 
multiplier  of  1 .  If  any  of  the  L/D  values  proved  to  be  valid,  the  pair  of  velocity  curves 
should  straddle  the  theoretical  velocity  curve  line.  Once  an  L/D  value  was  chosen,  the 
next  step  in  the  simulation  set-up  could  be  determined. 

The  second  step  in  the  experimentation  was  to  determine  the  desired  force  multiplier  for 
our  selected  L/D,  providing  an  experimental  velocity  curve  capturing  the  physics  of  the 
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problem.  Again  using  50  m/s  and  a  tube  radius  of  20A,  the  desired  force  multiplier  was 
experimentally  found.  The  body  force  multiplier  was  varied  from  ~0  to  1  for  the  selected 
L/D.  The  resulting  velocity  curves  were  graphed  against  the  theoretical  curve  and  a 
visual  comparison  was  made.  Once  the  force  multiplier  was  determined,  the  research 
effort  moved  into  the  next  phase  of  simulations. 

The  next  phase  of  simulations  moved  to  simulating  a  coaxial  flow  with  and  without 
acoustic  interference.  When  simulating  the  effects  of  acoustics  versus  the  absence  of 
acoustics,  a  separate  simulation  set-up  was  employed.  The  tube  radius  was  increased  to 
65A,  which  was  determined  to  be  the  minimum  radius  possible  based  on  Kn  number  for 
the  purposes  of  the  research.  This  radius  allows  the  molecular  dynamic  simulation  to 
adequately  represent  continuum  solutions  (Branam,  2005).  Tube  length  was  maintained 
at  five  times  the  tube  diameter.  The  chamber  width  was  also  set  at  five  times  the  tube 
diameter.  The  tube  wall,  chamber  wall  and  coaxial  wall  temperatures  were  set  at  120K, 
300K  and  300K  respectively.  The  tube  centerline  velocity  was  maintained  at  50  m/s 
while  the  coaxial  jet  velocity  was  varied  in  separate  simulations  at  0,  50  and  300  m/s. 
Later  sections  discuss  the  coaxial  jet  sizing  in  more  detail. 

Flow  Parameters 

Flow  in  a  molecular  dynamic  simulation  is  modeled  as  individual  interacting  molecules. 
The  flow  in  this  research  is  modeled  using  nitrogen  molecules.  Nitrogen  molecules  are 
simple  to  model  and  they  allow  for  the  consideration  of  a  rotating  molecule  in  the  flow 
field.  The  temperatures  maintained  inside  the  inner  tube,  outer  tube,  and  chamber  are 
120  K,  300  K,  and  300  K  respectively.  The  densities  maintained  inside  the  inner  tube, 
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outer  tube,  and  chamber  are  564.79  kg/m3,  45.087  kg/m3,  and  45.087  kg/m3  respectively. 
The  pressure  modeled  in  the  flow  is  4.0  MPa. 

Equations  of  State 

In  molecular  dynamics,  the  thermodynamic  properties  of  a  flow  are  not  inputs  but  rather 
the  results  of  the  simulation.  Due  to  the  output  of  molecular  dynamics,  turbulent  flow 
requires  no  special  considerations,  as  Navier-Stokes  models  do.  The  velocity  distribution 
becomes  an  important  factor  when  comparing  the  results  to  the  actual  physics  of  a 
problem.  The  velocity  distribution  of  the  particles  agrees  with  the  Maxwell/Boltzman 
distribution  at  equilibrium  (Branam,  2005:15). 

Molecular  dynamics  establishes  the  temperature,  density,  and  volume  for  a  desired 
pressure.  The  temperature  at  any  given  time  in  the  simulation  is  very  important  for  the 
understanding  the  results  of  the  simulation.  The  relationship  for  absolute  temperature  is  a 
direct  product  of  kinetic  theory  (Branam,  2005:22). 

Reynolds  and  Knudsen  Number 

Reynolds  number  and  Knudsen  number  ranges  play  a  large  part  in  getting  the  correct 
scaling  for  the  simulations  in  order  to  provide  results  comparable  to  experimental 
research.  The  small  size  of  the  simulations  requires  some  scaling,  while  flow  conditions 
along  with  molecular  dynamic  considerations  drive  other  variables. 

The  simulated  flow  begins  in  the  tube  as  a  laminar  flow.  As  the  flow  moves  into  the 
chamber,  a  turbulent  flow  field  is  expected  due  to  flow  instabilities.  Reynolds  numbers 
for  experimental  turbulent  tube  flows  in  these  conditions  are  usually  between  10,000  and 
100,000.  These  high  Re  are  the  calling  card  of  turbulent  flows.  The  limits  on  the  number 
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of  particles  molecular  dynamic  simulations  can  model  forces  Re  for  this  research  effort  to 
be  closer  to  1 .0.  The  low  Re  is  the  largest  indicator  of  uncertainty  when  modeling  flow 
fields  in  this  manner.  Regardless,  previous  work  has  shown  good  agreement  with 
experimental  results  (Branam,  2005). 

Knudsen  number  plays  a  larger  part  in  the  description  of  the  flow  as  well  as  the  acoustic 
interference.  For  a  slip  condition,  the  Kn  will  be  in  the  range  of  0.0 1  and  0.1.  The 
simulations  presently  are  operating  in  the  no-slip  regime  and  require  a  Kn  <  0.01. 

Coaxial  Jet  Sizing 

The  coaxial  jet  experimentation  was  to  be  compared  to  the  previous  dissertation  work  of 
Davis  (2006).  While  Davis’  work  was  qualitative,  this  research  will  be  both  a  quantitative 
and  qualitative  effort.  To  get  the  best  research  comparisons  for  the  coaxial  jet  flow,  the 
molecular  dynamic  model  needed  to  be  scaled  similarly  in  order  to  show  comparisons  to 
the  dimensions  of  the  physical  model  in  Davis’  work.  The  physical  model  appears  in 
Figure  5.  The  inner  jet  diameter  was  0.5  mm  while  the  outer  jet  diameter  was  2.42  mm. 
The  coaxial  jet  outer  to  inner  diameter  ratio  D0/Dj  was  4.84.  The  thickness  between  the 
inner  and  outer  core  was  give  as  0.54  mm. 
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Figure  5.  Close-up  image  of  the  coaxial  injector  tip  area,  the  dimensions  of  the  injector 
are:  Ri  =  0.25  mm,  R2  =  0.79  mm,  R3  =  1.21  mm  and  R4  =  1.59  mm.  (Davis  and 

Chehroudi,2006) 

Molecular  dynamic  simulation  is  limited  to  nanometers  as  opposed  to  mm.  The  coaxial 
tube  needs  to  be  scaled  while  maintaining  the  geometry  of  the  physical  model  as  much  as 
possible.  In  order  to  achieve  a  valid  molecular  simulation,  the  inner  tube  diameter 
needed  to  be  65A.  Using  this  initial  dimension  as  a  base,  the  ratios  of  the  physical  model 
determined  the  remaining  geometry  of  the  coaxial  jet.  The  resulting  geometry  gave  an 
inner  diameter  of  65 A,  an  outer  diameter  of  3 14. 6 A  and  140. 4 A  spacing  between  the 
inner  and  outer  jets.  This  geometry  maintained  the  D0/Dj  at  4.84.  Figure  6  provides  a 
pictorial  representation  of  the  molecular  dynamic  simulation  of  the  coaxial  jet. 
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Figure  6.  Pictorial  representation  of  the  coaxial  injector  tip  area,  the  dimensions  of  the 
injector  are:  Ri  =  65A,  R2  =  205. 4A,  R3  =  314.6A. 

After  running  multiple  simulations  using  the  new  geometry,  the  scaling  of  the  coaxial  jet 

needed  to  be  changed.  The  ratio  of  the  jet’s  outer  to  inner  diameter  caused  problems  as 

the  calculations  exceeded  the  available  resources  available  for  the  simulation.  The 

coaxial  jet  would  need  to  be  smaller  to  fit  the  capabilities  of  the  available  computing 

facilities.  The  ratio  of  outer  to  inner  diameter  D0/Di  for  the  physical  model  of  the  coaxial 

jet  started  at  4.84  but  needed  to  be  decreased  to  4.0.  D0/D,  was  the  only  ratio  change 

made  to  the  simulated  coaxial  jet.  The  inner  jet  diameter  was  maintained  at  65  A  and  the 

spacing  between  the  outer  and  inner  jet  was  maintained.  In  essence,  the  outer  jet  area 

was  decreased  to  fit  computing  capabilities.  Figure  7  provides  a  pictorial  representation 

of  the  new  simulated  coaxial  jet.  The  decrease  in  the  outer  jet  area  is  clearly  visible  as 

the  outer  jet  thickness  is  halved  from  109.2A  to  54. 6A. 
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Figure  7.  Pictorial  representation  of  the  coaxial  injector  tip  area,  the  dimensions  of  the 
injector  are:  Rj  =  65 A,  R2  =  205. 4A,  R3  =  260 A. 

While  the  coaxial  jet  geometry  and  inner  jet  velocity  will  not  change  throughout  the 

several  planned  simulations,  the  same  can  not  be  said  for  the  outer  jet  velocity.  While  the 

inner  jet  velocity  was  50  m/s,  the  outer  jet  velocity  will  be  varied  to  provide  a  wide  range 

for  Uo/Ui.  Due  to  the  present  simulation  run  times,  the  research  will  initially  look  at  three 

outer  jet  velocities.  The  outer  jet  velocities  will  be  0,  50  and  300  m/s  providing  a  U0/Ui 

of  0,  1,  and  6  respectively. 

Weber  number  (We)  is  a  dimensionless  ratio  of  momentum  force  to  surface  tension  force, 
which  helps  describe  shear  flow.  Weber  number  for  a  single  jet  is  expressed  in  equation 
43. 


We  = 


pU'd 


(43) 


cr 

When  calculating  We  for  a  coaxial  flow,  the  difference  between  inner  and  outer  jets  must 
be  accounted  for  and  the  equation  is  expanded  to 
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We  = 


(44) 


P(u:--u?)d 

a 

where  d  is  the  inner  tube  diameter  and  a  is  the  surface  tension.  The  simulation  does  not 
employ  a  two-phase  flow,  so  surface  tension  is  not  a  significant  factor.  Momentum  ratios 
and  momentum  flux  ratios  provide  similar  information  about  jets  without  taking  into 
account  surface  tension.  There  is  no  momentum  ratio  or  momentum  flux  of  outer  to  inner 
jet  (equations  45  and  46  respectively)  for  the  single  jet  case. 

M„=(p„U»2A)/(pV-A,)  (45) 

M  =(p,U>)l(p,U,2)  (46) 

Momentum  flux  ratio  is  a  key  quantity  for  near-field  geometry  of  a  coaxial  configuration 
(Villermaux,  1998:807).  Momentum  flux  ratio  highly  affects  jet  shapes  and  spray  angles. 
The  inner  potential  core  length  over  inner  diameter  ratio  shares  a  strong  relationship  with 
M.  For  single-phase  flow,  which  is  the  focus  of  this  research,  the  two  are  inversely 
related  with  inner  core  length  over  inner  diameter  ratio  ~  M  .  Due  to  this  inverse 
relationship,  very  low  momentum  ratios  (M<1)  cause  the  inner  jet  to  destabilize  very 
slowly,  increasing  the  inner  potential  core  length  to  inner  diameter  ratio  (Villermaux, 
1998:81 1).  The  momentum  ratios  for  U0/Ui  =  1  and  6  are  0.798  and  2.874,  respectively. 
While  the  jet  velocity  ratio  provides  information  about  the  flow,  momentum  ratio 
encompasses  the  flow  and  shape  of  the  jet  issuing  into  the  chamber. 

Chamber  Methodology 

After  the  flow  exits  the  tube,  it  is  issued  into  a  chamber  having  many  of  its  attributes 
defined  by  the  molecular  dynamic  simulation.  The  chamber  is  a  box  with  the  size 


33 


determined  by  the  user  with  fail-safe  mechanisms  correcting  the  sizing  if  it  is  out-of 
bounds  of  the  physical  requirements.  The  chamber  utilizes  periodic  boundary  conditions. 
For  data  analysis  purposes,  the  chamber  is  separated  into  eleven  segments.  Each  of  the 
segments  averages  the  properties  of  the  flow  field  within  that  segment  and  outputs  data 
accordingly.  The  data  profiles  are  displayed  radially  to  best  represent  the  axisymmetric 
nature  of  this  flow  problem. 

Wall  Methodology 

The  walls  of  the  tube  are  an  important  aspect  of  the  molecular  dynamic  simulation.  The 
molecular  dynamic  simulation  code  is  set-up  to  handle  multiple  wall  conditions.  Branam 
(2005:39-48)  has  done  work  to  compare  the  modeling  of  walls  as  a  system  of  particles. 
The  models  were  compared  to  a  diffusively  reflecting  wall  model.  The  diffusely 
reflecting  wall  represents  the  simplest  of  the  wall  models.  The  model  uses  a  Maxwell- 
Boltzmann  distribution  to  detennine  the  velocity  and  magnitude  of  the  particles  reflecting 
from  the  wall.  The  one  zone  wall  utilizes  springs  to  anchor  the  wall  particles  to  the 
lattice  site  and  includes  a  stochastic  zone.  “The  stochastic  zone  acts  like  a  thennal 
reservoir  to  dissipate  energy  transferred  to  the  wall  from  the  fluid  (Branam,  2005:41).” 
The  two-zone  wall  allows  the  wall  particles  to  vibrate  around  their  lattice  sites  while 
continuing  to  absorb  energy  from  the  fluid.  The  three-zone  wall  includes  a  reaction  zone 
directly  interacting  with  the  fluid  particles.  It  also  includes  a  stochastic  zone  and  a  layer 
of  stationary  particles.  This  research  relied  on  the  three-zone  wall  due  to  its 
demonstrated  ability  to  best  capture  the  physics  of  the  simulations.  The  three-zone  wall 
is  the  most  complex  of  the  three  wall  models  and  does  increase  computing  times.  Figure 
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8  shows  the  geometry  of  the  three-zone  wall.  The  three-zone  wall  uses  the  Lennard- 


Jones  potential  to  determine  the  intermolecular  forces  within  the  wall,  as  well  as  to 
govern  the  interactions  between  wall  and  flow  particles  (Branam,  2005:41-46). 
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Figure  8.  Three-Zone  Wall  Model  (Branam,  2005:46) 

Acoustics  Set-Up 

Acoustic  interference  is  simulated  in  two  manners  for  this  effort.  The  first  method  will 
simulate  the  wall  as  a  row  of  modeled  particles  applying  pressure  against  a  chamber  wall. 
The  second  method  will  model  the  acoustic  interference  using  a  force  model  representing 
discrete  wall  particles.  The  entire  molecular  dynamics  system  is  modeled  as  particles, 
including  walls  and  coaxial  structures.  Particle  modeling  was  the  inspiration  behind  the 
first  method  of  acoustic  interference.  Rows  of  particles  were  used  to  model  the 
interference  along  a  face  of  the  chamber.  The  molecules  were  given  amplitude  and 
frequency  in  order  to  mimic  the  desired  acoustic  set-up.  The  acoustic  interference  was 
modeled  for  an  SPL  of  90  dB  and  for  a  range  of  frequencies  between  10  and  12  kHz.  The 
frequency  scaling  required  to  affect  the  simulations  is  discussed  later  in  the  section. 

An  acoustic  frequency  was  selected  to  provide  comparable  results  to  previous 
experimental  research  efforts.  The  frequency  is  used  in  the  wave  equation  to  give  the 
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wave  displacement  amplitude  (“Decibel  Scales”,  undated).  The  pressure  amplitude  is 
then  calculated  and  compared  to  atmospheric  pressure  using  the  SPL  equation  producing 
a  decibel  level.  The  acoustic  molecules  were  modeled  using  SPL.  The  acoustic 
molecules  were  arranged  in  three  rows  and  placed  against  the  y-axis  of  the  chamber.  The 
molecules  could  then  affect  the  chamber  particles  by  moving  with  the  desired 
displacement  creating  a  force  against  the  chamber  particles. 

This  research  also  explored  a  mathematical  acoustic  interference  model.  If  the  acoustic 
pressure  necessary  to  provide  the  desired  displacement  could  be  decided  upon,  then  it 
would  be  possible  to  calculate  the  desired  force  from  the  definition  of  pressure. 

Similarly,  acoustic  pressure  can  be  added  to  the  force  calculations  based  on  the  definition 
of  pressure  felt  between  interacting  molecules.  The  pressure  is  calculated  from  the  virial 
work  term. 

Considering  a  virial  term  for  acoustic  pressure  analogous  to  the  term  used  to  define 
pressure,  an  acoustic  force  factor  (ffac0us )  results.  The  pressure  term  as  determined  from  a 
Lennard-Jones  potential  is 


Plj  ~ 


8g 

Vol 


f  _12 


12 


.6  A 


2  r 


.  ffuL 

6VoI 


(47) 


Substituting  an  acoustic  pressure  into  this  expression  and  solving  for  an  acoustic  force 
factor  yields 


ff  =- 

JJ  acous 


(6  -Vol) 


(48) 


r 

This  method  due  to  its  inverse  relationship  with  distance  between  interacting  molecules 
requires  the  placement  of  a  phantom  line  in  between  the  acoustic  particles  and  the 
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chamber  particles  to  capture  the  true  physics  of  the  problem  and  provide  a  computational 


limiter.  Once  an  acoustic  force  is  calculated,  it  needs  to  be  given  a  sinusoidal  component. 
The  sinusoidal  wave  representing  the  acoustic  force  then  becomes 

facous  =  ffacous  S™  (<Ot  ~  </>)  (49) 

The  operating  frequency  modeled  by  the  acoustic  interference  is  in  the  range  of  10.0  - 
12.0  kHz.  These  frequencies  do  not  affect  the  model  in  their  present  state  due  to  the  size 
of  the  simulations.  Frequency  scaling  is  needed  in  order  to  achieve  the  proper 
interference.  Frequencies  are  scaled  higher  to  affect  the  simulation.  A  reference  tube 
radius  of  1  cm  is  used  for  scaling  purposes.  The  operating  frequency  range  is  equated  to 
operating  on  the  reference  tube  radius.  The  frequency  is  then  scaled  higher  using 
equation  50. 

<x>  =  coo—  (50) 

r 


Parallel  Computing 

Due  to  the  computationally  intensive  calculations  required  for  molecular  dynamic 
simulations,  steps  are  taken  to  make  the  process  work  more  efficiently.  Running  the 
simulation  on  multiple  processors  leads  to  shortened  simulation  completion  times.  The 
computing  resources  consisted  of  a  cluster  of  networked  servers.  Sixty-four  nodes 
existed  offering  128  2.2  GHz  Opteron  248  CPUs.  Each  node  provides  4  GB  of  memory 
with  54  GB  of  local  work  partition  per  node.  The  nodes  were  also  capable  of  1  MB  of 
secondary  cache  per  processor.  The  parallelization  scheme  employed  Message  Passing 
Interface  (MPI)  in  the  code  to  take  advantage  of  the  available  computer  resources. 
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Increasing  the  number  of  processors  used  to  run  a  simulation  in  general  will  decrease  the 
wall  calculation  times  but  increases  the  overall  CPU  hours.  Branam  (2005)  found  as  the 
number  of  processes  used  increases  the  communication  time  also  increases  and  a 
breakeven  point  is  eventually  reached.  The  best  parallel  processing  compromise  for  the 
simulations  in  this  research  was  8  nodes  with  two  processes  per  node.  Branam  (2005) 
also  showed  the  utilization  of  the  velocity  Verlet  algorithm  and  RATTLE  constraints 
required  one  major  communication  per  time  step.  For  small  simulations  using  argon  and 
nitrogen  communication  times  equated  to  50%  of  the  simulation  run  time  on  16 
processes. 

The  code’s  parallelization  technique  used  particle  decomposition.  The  technique 
balances  calculations  between  processors  improving  efficiencies.  With  particle 
decomposition,  each  processor  is  responsible  for  the  calculations  of  a  set  of  particles 
rather  than  passing  particles  between  processors.  Each  processor  will  be  aware  of  all  the 
particles  but  will  only  be  responsible  for  the  force  calculations  of  their  specific  particles. 
This  method  is  a  fairly  efficient  technique  used  to  balance  force  and  neighbor 
calculations  but  it  reaches  an  upper  limit  near  two  million  particles  when  2.0  GB  per 
process  of  memory  is  available. 
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IV.  Results 


Chapter  Overview 

This  chapter  presents  the  results  of  the  research.  The  first  section  looks  at  the  initial  steps 
taken  to  match  the  theoretical  radial  velocity  profile  in  the  tube.  The  variation  of  the  tube 
length  over  tube  diameter  is  discussed  along  with  the  force  factor  variations  made  to 
match  the  theoretical  curve.  This  section  also  shows  the  actual  velocity  profiles  for  the 
tube  and  coax  produced  during  the  simulations.  The  simulated  velocity  curves  are 
compared  with  the  theoretical  curves  and  conclusions  are  drawn.  The  next  section  delves 
into  the  acoustic  models  used  in  this  research.  The  two  acoustic  interference  models  are 
compared  to  cases  with  the  acoustics  turned  off  in  order  to  validate  each  of  the  models. 
Next,  mixing  along  the  centerline  of  the  chamber  is  investigated.  The  effects  of  acoustic 
interference  on  mixing  are  discussed  for  all  three  cases.  The  effects  of  different  velocity 
ratios  on  the  mixing  are  also  investigated.  The  radial  density  profiles  for  each  case  are 
examined  in  the  final  section.  The  density  profile  provides  valuable  information  about 
the  shape  of  the  jets  issuing  into  the  chamber. 

Matching  Theoretical  Velocity  Profile 

The  initial  experimentation  used  to  detennine  the  minimum  tube  length  necessary  to  get  a 
valid  flow  profile  was  completed  using  an  earlier  and  simpler  molecular  dynamic 
simulation  tool.  In  order  to  detennine  the  valid  velocity  profile  matching  the  theoretical 
flow  profile,  the  simulation  ran  multiple  times  varying  both  the  ratio  of  tube  length  to 
tube  diameter,  L/D  and  the  force  factor  used  to  create  the  Poiseuille  flow  conditions 
throughout  the  tube  and  chamber.  The  injector  simulation  was  set  to  an  initial  velocity  of 
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50  m/s.  For  a  fixed  diameter  size  of  20  A,  injection  was  simulated  for  L/D  values  of  1,  5 
and  10.  Each  L/D  value  was  also  simulated  for  a  force  factor  of  approximately  0  and  1. 
The  ideal  case  would  find  the  theoretical  flow  profile  line  straddled  by  the  experimental 
flow  profiles  for  a  force  factor  of  ~0  and  1.  Figures  9,  10  and  1 1  show  the  flow  profiles 
of  L/D  =1,10  and  5  respectively  as  they  are  compared  to  the  same  theoretical  velocity 
profile.  Velocity  profile  charts  do  not  represent  velocity  as  it  flows  through  the  tube. 

The  x-axis,  r/rtube,  of  the  velocity  profile  is  representative  of  the  change  from  0.0,  which 
represents  the  centerline  of  the  inner  tube  to  1.0,  which  represents  the  wall  of  the  inner 
tube. 

Analytical  flow  profiles  for  L/D  =  1  at  force  factors  of  ~0  and  1  did  not  surround  the 
theoretical  flow  profile  and  were  discarded  as  a  potential  flow  geometry.  Llow  profiles 
for  L/D  =  5  and  10  both  fit  the  criteria  for  possibly  matching  the  theoretical  flow  profile. 
Lor  these  two  geometries,  the  theoretical  flow  profile  sits  between  the  force  factor 
profiles.  Since  both  are  feasible  geometries,  L/D  =  5  was  selected  to  minimize 
computational  run  time. 
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Figure  9.  Radial  Velocity  Profile,  rtabe  =  20A,  L/D  =  1 
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Figure  10.  Radial  Velocity  Profile,  rtube  =  20A,  L/D  =  10 
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A  force  factor  of  ~0  provided  a  flow  profde  below  the  desired  initial  speed  of  50  m/s 
while  a  force  factor  of  1  produced  an  initial  flow  velocity  almost  twice  the  theoretical 
value.  The  force  factor  was  then  adjusted  by  0.05  from  0.05  to  1  until  the  analytical 
profile  matched  closely  the  theoretical  profile.  The  curve  alignment  occurred  at  a  force 
factor  of  0.4  (Figure  1 1).  The  final  geometry  chosen  was  L/D  =  5  with  a  force  factor  of 
0.4.  This  geometry  and  the  force  factor  were  maintained  throughout  the  rest  of  the 
simulations. 

Velocity  Profile 

Figure  12  shows  the  velocity  profile  for  the  single  jet  flow.  The  velocity  profile  error  bars 
represent  a  spatial  resolution  that  captures  the  values  and  trends  of  the  data.  The  velocity 
curve  for  acoustic  free  and  acoustic  interference  models  show  acceleration  in  the  tube 
after  the  simulation  is  first  begun.  This  same  trend  appears  later  in  the  other  simulation 
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cases.  The  transience  caused  by  starting  the  jet  to  issue  into  the  chamber  causes  the 
acceleration  of  the  tube  flow.  The  flow  of  molecules  appears  to  rush  into  the  chamber 
but  as  the  simulation  moves  towards  the  tube  walls,  the  velocity  decreases  until  it 
converges  closely  to  the  expected  value.  The  single  jet  was  previously  studied  and 
matching  the  theoretical  curve  has  been  accomplished  but  more  may  need  to  be  done  to 
ensure  the  coaxial  flow  field  can  be  controlled  similarly. 


Figure  12.  Radial  Velocity  Profde  for  U0  =  0  m/s,  U;  =  50  m/s 
Figures  13  and  14  show  the  trends  of  both  the  acoustic  free  and  acoustic  interference 
models  for  U0/Ui  =  1  and  6  respectively.  The  addition  of  the  coaxial  component  to  the 
flow  creates  some  interesting  results.  For  U0/Uj  =  1  and  6,  similar  results  can  be  seen. 
The  theoretical  curve  for  the  inner  tube  remains  unchanged  but  the  acoustic  free  velocity 
curves  for  each  show  increased  initial  velocity.  The  two  acoustic  free  velocity  curves  for 
U0/Ui  =  1  and  6  show  an  increased  initial  velocity  in  the  tube.  The  acoustic  interference 
curve  shows  acceleration  in  both  coaxial  cases.  In  Figure  13,  the  acoustic  free  model  for 
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U0/Ui  =  1  shows  smaller  levels  of  acceleration  when  compared  to  the  other  cases.  The 
probable  reason  is  the  level  of  data  collected  for  this  case  is  small  when  compared  to 
other  cases.  Acceleration  is  to  be  expected  in  the  acoustic  free  case  for  U0/Ui  =  1  once 
more  data  is  collected.  The  acoustic  free  and  acoustic  interference  models  for  the  U0/Ui  = 
6  case  have  a  very  similar  centerline  velocity,  leading  to  the  conclusion  the  acoustic 
interference  is  not  the  leading  effect  in  starting  the  tube  flow.  The  addition  of  the  coaxial 
jet  affects  the  pressure  for  flow  issuing  out  of  the  tube  centerline. 


Figure  13.  Radial  Velocity  Profde  for  U0  =  50  m/s,  U;  =  50  m/s 
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Figure  14.  Radial  Velocity  Profile  for  U0=  300  m/s,  If  =  50  m/s 


A  theoretical  velocity  profile  for  the  outer  tube  was  calculated  for  comparison  to  the 
simulated  coaxial  profiles.  The  velocity  profile  charts  for  the  outer  tube  differ  from  the 
velocity  profile  charts  for  the  inner  tube.  The  y-axis  still  represents  velocity  while  the  x- 
axis  still  represents  a  radius  ratio.  The  x-axis,  r/(r0-ri),  represents  the  inner  wall  of  the 
outer  tube  at  0.0  to  the  outer  wall  of  the  outer  tube  at  1.0.  Figures  15  and  16  show  the 
theoretical  coaxial  velocity  curves  for  U0  =  50  m/s  and  300  m/s  respectively.  The 
theoretical  coaxial  velocity  curve  for  both  cases  is  plotted  across  the  entire  coaxial  radius. 
The  curve  across  the  coaxial  radius  takes  a  parabolic  shape  peaking  close  to  the 
centerline.  The  simulated  coaxial  curves  do  not  in  any  way  follow  the  trend  of  the 
theoretical  curve.  The  curves  for  the  acoustic  free  and  acoustic  interference  models 
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follow  similar  trends  especially  in  the  U0/Ui  =  6  case.  These  results  need  further 


investigation. 


r/(r0-r;) 

Figure  15.  Coaxial  Velocity  Profde  for  U0  =  50  m/s,  Ui  =  50  m/s 
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The  simulated  coaxial  curve  for  U0  =  50  m/s  does  initially  begin  at  50  m/s  but  shows  a 
large  drop-off  in  velocity  before  accelerating  again  only  to  decelerate  further  from  the 
centerline  with  the  curve  ending  in  an  acceleration.  The  shape  of  the  curve  resembles  a 
decaying  sinusoidal  wave.  The  coaxial  velocity  profde  for  U0  =  300  m/s  has  a  somewhat 
similar  profile  but  initially  the  flow  in  this  case  is  moving  at  a  negative  rate  and  begins  to 
accelerate  fluctuating  until  it  balances  towards  the  outer  coaxial  wall.  As  discussed 
earlier,  the  size  of  the  coaxial  component  had  to  be  decreased  to  stay  within  the 
computational  limits  of  the  computer  network.  The  shape  of  the  velocity  profiles  is  most 
likely  due  to  small  gap  width  in  the  coaxial  tube.  As  stated  earlier,  Kn  for  continuum 
flow  should  be  less  than  0.01.  Kn  reached  0.29  in  the  coaxial  component  which  is  well 
outside  the  limits  of  continuum  flow.  The  decreased  size  of  the  coaxial  jet  may  have 
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affected  the  simulation.  The  smaller  size  does  not  seem  to  be  able  to  accurately  capture 
the  physics  of  the  problem. 

Acoustic  Model  Validation 

Mass  fraction  charts  represent  the  mixing  of  the  atoms  from  the  tube  and  coax  with  the 
chamber  particles.  Lower  mass  fractions  represent  higher  mixing  in  the  chamber.  The 
mass  fraction  charts  show  mass  fractions  ranging  from  0.0  to  1.0  on  the  y-axis.  The  x- 
axis,  z/d,  is  the  non-dimensionalized  chamber  length.  It  was  stated  earlier  that  the 
chamber  is  divided  into  1 1  sections.  The  non-dimensionalized  values  of  the  x-axis 
represent  the  1 1  chamber  sections.  Figure  17  shows  the  earliest  captured  mixing  at  U0/Ui 
=  6  for  all  three  models.  Even  at  this  early  stage,  there  is  much  that  can  be  drawn  from 
the  chart.  The  acoustic  free  model  and  discrete  particle  model  produce  very  similar 
curves.  The  force  model  curve  on  the  other  hand  shows  increased  mixing  in  the  chamber. 
Observing  the  flow  after  half  of  a  nanosecond  provides  similar  results.  The  acoustic  free 
model  and  the  discrete  particle  model  appear  to  lie  on  top  of  each  other.  The  force  model 
again  shows  increased  mixing  when  compared  to  the  acoustic  free  model  as  the  flow 
develops.  As  the  flow  develops  in  time,  the  force  model  mixing  seems  to  decrease  and 
converge  to  the  level  of  mixing  in  the  acoustic  free  model.  The  force  model  then 
separates  and  produces  lower  mixing  fractions  than  the  acoustic  free  case.  Even  though 
Figures  17  and  18  represent  flows  that  have  not  reached  a  steady  state,  from  these  figures, 
it  is  the  discrete  particle  model  of  acoustic  interference  is  not  affecting  the  chamber  as 
intended.  The  probable  reason  for  the  ineffectiveness  of  the  model  is  the  need  for  a 
larger  system  and  more  acoustic  particles.  The  acoustic  particles  in  the  discrete  particle 
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model  presently  operate  at  a  90dB  SPL.  At  90  dB,  the  Lennard- Jones  forces  for  the 
acoustic  atoms  are  negligible  after  9.89  A.  The  mean  free  path  for  our  chamber  is 
approximately  16  A.  Due  to  the  size  scale  of  the  simulation,  simply  increasing  the 
decibel  level  to  100  dB  creates  a  20  A  acoustic  displacement,  too  big  to  affect  any 
chamber  particles.  The  present  decibel  level  is  limited  to  a  regime  unable  to  affect  our 
chamber  at  the  present  sizes. 


Figure  17.  Centerline  Mixing  Comparison  at  0.25  ns,  U0  =300  m/s,  Ui  =  50  m/s 
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Figure  18.  Centerline  Mixing  Comparison  at  0.75  ns,  U0  =300  m/s,  Ui  =  50  m/s 
Figures  19  and  20  show  the  acoustic  free  and  force  model  cases  for  a  coaxial  velocity  of 
50  m/s.  The  acoustic  discrete  particle  model  is  omitted  from  the  graph  since  it  does  not 
capturing  the  physics  of  the  problem.  Figures  19  and  20  show  similar  mass  fraction 
trends  validating  the  force  model  of  acoustic  interference.  Increased  mixing  is  evident  in 
the  force  model  case  in  both  charts.  The  acoustic  free  model  does  not  exhibit  mixing  on 
the  level  seen  for  the  acoustic  interference  model. 
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Figure  19.  Centerline  Mixing  Comparison  at  0.25  ns,  U0  =50  m/s,  U,=  50  m/s 


Figure  20.  Centerline  Mixing  Comparison  at  1.0  ns,  U0  =50  m/s,  U,  =  50  m/s 
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Centerline  Flow  Development 

The  development  of  the  flow  in  the  chamber  is  an  important  aspect  in  verifying  the 
simulation.  Looking  at  the  six  cases  and  following  the  development  of  the  flow  through 
centerline  mass  fraction  data  provides  a  good  look  into  how  the  flow  is  mixing 
throughout  the  chamber.  The  mass  fraction  is  plotted  against  the  non-dimensionalized 
axial  chamber  dimension.  The  mass  fraction  is  the  ratio  of  mass  injected  into  the 
chamber  over  the  mass  at  the  specified  location  in  the  flow.  Again,  the  mass  fraction 
gives  insight  into  the  level  of  mixing  in  the  chamber. 

The  first  two  cases  to  be  discussed  are  U0  =  0  m/s  with  and  without  acoustic  interference. 
Figure  21  shows  the  development  of  the  flow  at  the  centerline  for  U0  =  0  m/s  for  the 
acoustic  free  case.  The  flow  develops  as  time  progresses,  steadily  up  to  about  1.0  ns.  As 
time  continues  past  1 .0  ns,  the  flow  mixing  increases  as  the  flow  moves  longitudinally 
towards  the  center  of  the  chamber.  The  variations  seen  in  mass  fraction  ratios  results 
from  the  initial  acceleration  of  the  jet  into  the  chamber.  Mass  fraction  data  at  2.0  and  3.0 
ns  shows  decreased  mixing  when  compared  to  the  1.0  ns  curve  but  the  data  also  shows 
the  flow  moving  close  to  a  steady  state.  As  the  flow  moves  further  into  the  chamber,  the 
flow  exhibits  a  slight  oscillatory  nature.  This  oscillatory  motion  is  apparent  between  2.0 
ns  and  3.0  ns. 

Figure  22  shows  the  effects  of  acoustic  interference  on  the  flow  when  U0  =  0  m/s.  The 
flow  acts  similar  to  the  acoustic  free  flow.  The  mass  fraction  increases  with  time  and 
reaches  a  point  at  approximately  2.0  ns  where  the  flow  mixing  begins  to  converge. 
Looking  at  2.0  and  3.0  ns  curves  it  appears  the  flow  in  the  chamber  near  the  injector  has 
reached  a  steady  state  condition.  Figure  23  shows  a  running  comparison  of  the  acoustic 
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free  case  with  the  acoustic  interference  case.  At  the  early  stages  of  the  simulation,  very 
little  difference  is  observed  between  the  acoustic  free  and  acoustic  interference  cases.  As 
the  simulation  time  increases  and  the  flow  moves  further  into  the  chamber,  the  acoustic 
free  model  for  U0=  0  m/s  exhibits  more  mixing  than  the  acoustic  interference  model.  The 
acoustic  free  case  mixes  the  injected  mass  much  better  throughout  the  chamber  and 
shows  much  higher  mass  dissipation.  The  jet  reaches  as  far  as  z/d  ~  10  as  evidenced  by 
the  presence  of  mass  injected  into  the  chamber.  The  acoustic  free  model  also  seems  to 
move  closer  to  a  steady  state  value  at  3.0  ns.  The  mass  fraction  curves  at  2.0  and  3.0  ns 
show  little  difference  while  the  acoustic  interference  case  needs  longer  simulation  times 
to  reach  a  steady  condition. 


Figure  21.  Mixing  Along  Centerline  ,  Uo=0  m/s,  Acoustic  Free  Model 
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Figure  22.  Mixing  Along  Centerline  ,  Uo=0  m/s,  Acoustic  Interference  Model 


Figure  23.  Mixing  Along  Centerline,  Comparison  of  Acoustic  Free  Model  vs.  Acoustic 

Interference  Model  for  Uo=0  m/s 
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Figure  24  shows  the  centerline  mixing  for  the  U0  =  50  m/s  case  at  1.0  ns  for  the  acoustic 
free  case.  Presently,  1 .0  ns  is  the  extent  of  the  data  gathered  for  this  case.  This  chart 
does  not  provide  much  data,  but  the  data  it  does  present  is  compared  to  other  cases  later 
in  the  section.  Figure  25,  which  relates  to  the  acoustic  interference  case  for  U0  =  50  m/s 
shows  similar  mixing  as  the  flow  develops  over  time.  The  mixing  is  again  higher  at  1 .0 
ns  and  decreases  as  it  reaches  2.0  ns.  This  effect  is  attributed  to  the  initial  acceleration 
into  the  chamber.  The  flow  develops  as  simulation  time  increases  but  the  limited  data 
shows  the  flow  mixes  better  in  the  acoustic  interference  case  than  it  does  in  the  acoustic 
free  case.  The  shape  of  the  mixing  curves  in  the  acoustic  free  case  through  the  chamber 
is  similar  to  the  acoustic  interference  but  the  mixing  is  definitely  higher  for  the  acoustic 
interference  model.  Figure  26  shows  this  comparison  more  clearly.  For  the  earlier 
simulation  times,  it  is  clear  the  acoustic  free  mixing  is  lower. 


Figure  24.  Mixing  Along  Centerline,  Uo=50  m/s,  Acoustic  Free  Model 
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Figure  25.  Mixing  Along  Centerline,  Uo=50  m/s,  Acoustic  Interference  Model 


Figure  26.  Mixing  Along  Centerline,  Comparison  of  Acoustic  Free  Model  vs.  Acoustic 

Interference  Model  for  Uo=50  m/s 
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Figures  27  to  29  relate  to  U0  =  300  m/s.  The  development  of  this  flow  closely  resembles 
the  development  of  the  flow  for  U0  =  50  m/s.  For  the  acoustic  free  case,  the  flow 
develops  as  expected  with  mass  fraction  increasing  throughout  the  chamber  as  the 
simulation  time  increases.  The  acoustic  interference  case  faces  similar  development.  As 
the  simulation  time  increases,  mass  fraction  of  the  flow  increases.  The  flow  begins  to 
exhibit  some  signs  of  possibly  reaching  a  steady  state  flow  as  the  simulation  time  moves 
into  the  2.0  to  3.0  ns  region.  The  mixing  curves  begin  to  level  off  and  look  very  similar. 
The  difference  between  the  curves  at  this  time  may  be  due  to  the  acoustic  pressure  wave 
being  applied  to  the  flow.  This  pressure  could  cause  a  continued  fluctuation  in  mixing  as 
the  flow  simulation  time  increases.  Again,  as  in  the  U0  =  50  m/s  case,  the  mass  fraction 
in  the  acoustic  free  simulation  is  greater  than  seen  in  the  acoustic  interference  case. 
Similar  to  what  was  observed  for  U0  =  50  m/s  cases,  as  the  simulation  times  increase  the 
cases  with  and  without  acoustic  interference  show  some  correlation  in  mixing  early  in  the 
chamber  but  seem  to  diverge  as  the  flow  moves  further  into  the  chamber,  with  the  higher 
levels  of  mixing  in  the  acoustic  interference  model. 
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Figure  27.  Mixing  Along  Centerline,  Uo=300  m/s,  Acoustic  Free  Model 


Figure  28.  Mixing  Along  Centerline,  Uo=300  m/s,  Acoustic  Interference  Model 
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Figure  29.  Mixing  Along  Centerline,  Comparison  of  Acoustic  Free  Model  vs.  Acoustic 

Interference  Model  for  Uo=300  ms 

The  following  figures  show  comparisons  of  the  different  outer  jet  speeds  and  the  effects 
of  the  speeds  on  the  mixing  in  the  flow.  Figure  30  shows  all  three  outer  jet  velocity 
cases  at  1.0  ns  for  the  acoustic  free  model.  The  first  observation  made  is  the  U0  =  50  m/s 
and  U0  =  300  m/s  flows  seem  to  exhibit  very  similar  development  at  this  time  step  as  it 
will  be  seen  throughout  the  entire  simulation.  The  case  of  U0  =  0  m/s  exhibits  less 
injected  mass  in  the  acoustic  free  model.  Looking  at  Figure  31,  similar  observations  are 
expected  once  more  data  is  collected.  The  U0  =  50  and  300  m/s  cases  should  exhibit  very 
close  mass  fractions.  The  additional  mass  being  injected  by  the  coaxial  flow  as  well  as 
the  reduced  shear  forces  on  the  jet  should  account  for  the  higher  mass  fractions  seen  in 
the  U0/Ui  =  1  and  6  cases. 
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Looking  at  the  acoustic  interference  model,  some  interesting  observations  are  made.  For 
the  flow  at  1.0  ns  seen  in  Figure  33,  U0/Ui  =  0  and  6  cases  appear  to  exhibit  the  same 
mixing  in  the  chamber.  The  U0/Uj  =  1  case  exhibits  better  mixing  than  the  other  two 
cases.  At  2.0  ns  (Figure  34),  all  cases  begin  to  exhibit  different  levels  of  mixing 
throughout  the  chamber.  Figure  35  shows  as  simulation  time  increases  the  mixing  in  the 
chamber  for  the  U0  =  0  m/s  and  U0  =  300  m/s  cases  diverges  even  more.  The  mass 
fraction  throughout  the  chamber  at  3.0  ns  is  higher  for  U0  =  0  m/s  than  it  is  for  U0  =  300 
m/s,  as  one  would  expect,  due  to  the  higher  shearing  force  the  coaxial  flow  places  on  the 
core  jet  flow. 


Figure  30.  Mixing  Along  Centerline,  Comparison  at  1.0  ns  for  Acoustic  Free  Model 
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Figure  3 1 .  Mixing  Along  Centerline,  Comparison  at  2.0  ns  for  Acoustic  Free  Model 


Figure  32.  Mixing  Along  Centerline,  Comparison  at  3.0  ns  for  Acoustic  Free  Model 
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Figure  33.  Mixing  Along  Centerline,  Comparison  at  1.0  ns  for  Acoustic  Interference 

Model 
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34.  Mixing  Along  Centerline,  Comparison  at  2.0  ns  for  Acoustic  Interference 

Model 


Figure  35.  Mixing  Along  Centerline,  Comparison  at  3.0  ns  for  Acoustic  Interference 

Model 
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Radial  Density  Profile 

The  radial  density  profile  data  provides  some  insight  into  the  shape  of  the  jet  and  its 
issuance  into  the  chamber.  This  data  is  similar  to  that  of  a  jet  shadowgraph.  The  first 
simulation  discussed  is  the  single  jet  with  no  coaxial  component.  Again,  the  inner  tube 
velocity  for  all  cases  is  set  at  50  m/s.  Figure  37  shows  the  acoustic  free  model  of  the  jet 
soon  after  it  begins  issuing  into  the  chamber.  The  x-axis,  labeled  r/d,  is  a  non- 
dimensionalized  variable  equal  to  the  radial  chamber  position  divided  by  the  tube 
diameter.  The  y-axis  represents  a  density  ratio  calculated  as 

Density  Ratio  =  — — —  (51) 

Pt~  Pc 

The  t  and  c  subscripts  refer  to  tube  and  chamber,  respectively. 

At  z/d  =  0,  jet  exit  plane,  the  first  sign  of  tube  particles  issuing  into  the  chamber  can  be 
seen.  The  density  ratio  is  expected  to  be  close  to  1 .0  but  does  not  reach  this  level  (Figure 
36).  The  problem  is  most  likely  due  to  an  error  in  the  code.  Particles  are  not  being  added 
to  the  flow  at  the  expected  levels,  which  is  decreasing  the  density  ratio  seen  in  the 
chamber.  In  order  to  examine  the  radial  density  profiles  the  data  needs  to  be  scaled.  All 
radial  density  curves  are  scaled  to  the  first  point  near  the  exit  plane.  The  scaling  places 
the  first  data  point  at  1.0,  and  scales  all  other  points  in  accordance.  Figure  37  shows  the 
scaled  data  that  appeared  in  Figure  36.  The  first  curve  reaches  unity  as  initially  expected. 
Further  into  the  chamber  at  z/d  =  1,  the  presence  of  particles  can  be  seen  but  not  at  the 
levels  of  the  density  seen  at  the  jet  exit.  Moving  further  into  the  tube,  there  is  a  negligible 
density  change.  The  density  ratio  is  also  largest  at  the  center  of  the  jet  and  decreases  as  it 
moves  away  from  the  center  of  the  chamber  as  expected.  Figure  38  shows  the  acoustic 
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free  single  jet  flow  2.0  ns  into  the  flow  and  1.0  ns  after  the  events  in  Figure  37.  The  jet  is 
issuing  further  into  the  chamber  as  the  density  ratio  further  into  the  chamber  has 
increased  considerably.  The  expected  jet  growth  is  clearly  apparent  as  time  progresses. 
The  density  ratio  steadily  decreases  and  becomes  negligible  at  z/d  =  4.  Due  to  the 
scaling,  the  density  ratio  at  all  sections  of  the  chamber  can  be  compared  in  a  similar 
manner.  The  jet  is  changing  significantly  up  to  2.0  ns.  The  difference  in  density  ratios 
between  2.0  and  3.0  ns  near  the  exit  are  quite  small.  Greater  differences  are  seen  in 
density  ratio  between  z/d  =  2  and  z/d  =  4.  Though  changes  are  taking  place,  the  density 
ratios  seem  to  be  very  comparable  early  in  the  jet  and  it  leads  to  the  conclusion  the  jet 
ismoving  towards  steady  state. 


Figure  36.  Radial  Density  Profile  at  1.0  ns  without  Scaling,  U0  =  0  m/s,  Acoustic  Free 

Model 
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Figure  37.  Radial  Density  Profile  at  1.0  ns,  U0  =  0  m/s,  Acoustic  Free  Model 
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Figure  39.  Radial  Density  Profile  at  3.0  ns,  U0  =  0  m/s,  Acoustic  Free  Model 
The  time  development  of  the  density  ratio  follows  the  results  in  the  velocity  profiles  and 
mass  fraction  charts.  Figures  40,  41,  and  42  show  the  time  progression  of  density  at  z/d  = 
0,  1,  and  2  respectively.  The  time  progression  shows  an  initial  large  influx  of  particles 
into  the  chamber  for  the  first  1.0  ns,  with  a  decline  to  a  steady  state  in  the  times 
thereafter.  The  mass  fraction  time  profile  is  similar,  with  mixing  ramping  to  large  levels 
in  the  first  1.0  ns  and  then  decreasing  before  leveling  off  around  3.0  ns.  There  is  also 
some  correlation  with  the  velocity  profile.  The  flow  seems  to  have  some  initial 
acceleration  in  the  tube  but  the  flow  begins  to  decelerate  until  it  reaches  a  velocity  close 
to  the  theoretical  velocity. 
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igure  40.  Radial  Density  Profile  at  z/d  =  0,  U0=  0  m/s,  Acoustic  Free  Model 


Figure  41.  Radial  Density  Profile  at  z/d  =  1,  U0  =  0  m/s,  Acoustic  Free  Model 
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Figure  42.  Radial  Density  Profile  at  z/d  =  2,  U0  =  0  m/s,  Acoustic  Free  Model 
Looking  at  the  acoustic  interference  case  for  U0  =  0  m/s,  similar  trends  can  be  seen. 
Figure  43  shows  the  density  progression  at  1.0  ns.  The  acoustic  interference  model  does 
have  higher  density  ratios  than  the  acoustic  free  model  further  into  the  chamber.  Figures 
43,  44,  and  45  represent  1.0,  2.0  and  3.0  ns  after  the  start  of  the  flow  simulation.  The 
radial  density  for  the  acoustic  interference  model  follows  the  same  time  progression  as 
seen  in  the  acoustic  free  model.  The  biggest  observable  difference  is  the  higher  density 
ratios  at  each  section  of  the  chamber  exhibited  with  the  acoustic  interference  model. 
Figures  46  and  47  represent  time  progression  at  z/d  =  0  and  z/d  =  3  in  the  chamber.  The 
progression  follows  closely  to  the  acoustic  free  model.  Again,  the  density  ratio  increases 
as  time  grows  to  approximately  2.0  ns.  At  this  point,  the  density  ratio  begins  to  decrease 
until  it  begins  to  approach  a  steady  state  flow  at  3.0  ns. 
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Figure  43.  Radial  Density  Profile  at  1.0  ns,  U0  =  0  m/s,  Acoustic  Interference  Model 
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Figure  45.  Radial  Density  Profile  at  3.0  ns,  U0  =  0  m/s,  Acoustic  Interference  Model 


Figure  46.  Radial  Density  Profile  at  z/d=0,  U0  =  0  m/s,  Acoustic  Interference  Model 
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Figure  47.  Radial  Density  Profile  at  z/d=3,  U0  =  0  m/s,  Acoustic  Interference  Model 
Examining  the  affects  of  the  acoustic  interference  on  the  shape  of  the  jet  requires  side-by- 
side  comparisons  of  each  model  to  determine  what  part  acoustic  interference  plays  in  the 
jet  shape.  Figures  48  to  50  show  head-to-head  comparisons  of  the  single  jet  flow  at  1.0, 
2.0,  and  3.0  ns  respectively.  At  each  period,  the  acoustic  interference  model  shows 
higher  density  ratios.  Looking  at  the  work  of  Davis  (2006:94),  close  comparisons  can  be 
visually  drawn.  Figure  51  shows  experimental  work  studying  the  effects  of  acoustic 
interference  on  a  coaxial  jet.  Looking  at  the  first  column  in  Figure  51,  the  jet  in  frame  ‘a’ 
has  no  acoustic  interference  and  maintains  the  shape  expected.  Frame  ‘e’  shows  the  same 
jet  with  acoustic  interference  added.  The  first  observation  is  the  jet  has  become  radially 
thicker.  The  jet  also  has  become  more  unstable  and  does  not  hold  its  shape  as  well  as  the 
acoustic  free  jet.  The  experimental  results  suggest  wider  density  profdes  at  lower  values 
when  acoustic  interference  is  present.  The  acoustic  interference  simulations  are  showing 
higher  density  ratios  in  the  radial  axis.  The  density  ratios  for  the  acoustic  interference 
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simulations  also  show  higher  density  ratios  than  the  acoustic  free  model  further  into  the 
chamber.  These  results  do  not  seem  to  match  the  experimental  results  seen  but  are  most 
likely  explainable  by  the  jet  shape  and  discrete  particle  models  not  operating  within  the 
physics  of  the  problem. 


Figure  48.  Radial  Density  Profile  Comparison,  U0  =  0  m/s,  1.0  ns,  Acoustic  Free  (AF) 

Model  and  Acoustic  Interference  Model 
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Figure  49.  Radial  Density  Profile  Comparison,  U0  =  0  m/s,  2.0  ns,  Acoustic  Free  (AF) 

Model  and  Acoustic  Interference  Model 


r/d 

Figure  50.  Radial  Density  Profile  Comparison,  U0  =  0  m/s,  3.0  ns,  Acoustic  Free  (AF) 

Model  and  Acoustic  Interference  Model 
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Figure  5 1 .  Shadowgraph  images  of  coaxial  jet  at  subcritical  chamber  pressure  (~1 .5MPa) 
and  at  the  high  outer-jet  temperature  (~  190K)  corresponding  to  cases  1  -  4  (U0/Ui  =  ~  1, 
~2.5,  ~5,  and  —6.5).  The  acoustic  driver  is  off  for  images  in  the  top  row  and  on  for  the 
bottom  row  at  ~3  KHz.  The  velocity  ratio  is  the  same  for  each  column  and  increases 

from  left  to  right.  (Davis,  2006:94) 
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Figure  52  shows  the  density  profiles  for  the  U0  =  50  m/s  with  acoustics  turned  off.  The 
development  matches  closely  what  is  seen  in  the  previous  case.  Again,  the  data  collected 
are  very  preliminary  for  this  case  and  speaking  to  the  development  of  the  acoustic  free 
case  is  not  feasible.  Some  comparisons  can  be  made  between  the  acoustic  free  and 
acoustic  interference  cases  early  in  the  flow.  These  results  follow  later  in  this  section. 


Figure  52.  Radial  Density  Profile,  1.0  ns,  U0  =  50  m/s,  Acoustic  Free  Model 


Figures  53  and  54  show  density  ratio  development  of  the  acoustic  interference  case  for 
U0  =  50  m/s.  The  density  ratio  is  developing  similar  to  what  was  seen  in  the  U0  =  0  m/s 
case.  The  density  ratios  are  increasing  up  to  2.0  ns  and  the  flow  is  moving  further  into 
the  chamber.  Figure  53  shows  a  data  point  that  doesn’t  follow  a  smooth  trend  for  the  z/d 
=  0  line.  The  data  point,  at  a  value  of  r/d  =  0.5,  is  most  likely  due  to  a  data  collection 
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error  since  the  point  does  not  exist  in  similar  charts  at  different  velocity  ratios.  More  data 


is  necessary  to  assure  that  the  flow  is  moving  towards  a  steady  state. 


Figure  53.  Radial  Density  Profde,  1.0  ns,  U0  =  50  m/s,  Acoustic  Interference  Model 
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Figure  54.  Radial  Density  Profile,  2.0  ns,  U0  =  50  m/s,  Acoustic  Interference  Model 
Figures  55,  56,  and  57  show  development  at  z/d  =  0,  1,  and  2  respectively.  Each  chart 
shows  higher  density  ratios  at  2.0  ns  showing  a  growing  jet  up  to  this  time.  Figure  58 
shows  the  comparisons  of  the  acoustic  free  and  acoustic  interference  cases  for  the  U0  = 

50  m/s  case  at  1.0  ns.  These  results  match  what  is  seen  in  the  mass  fraction  charts.  This 
case  shows  higher  density  ratios  for  the  acoustic  free  case  which  is  a  departure  from  what 
is  seen  in  the  U0  =  0  m/s  case. 
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Figure  55.  Radial  Density  Profile,  z/d=0,  U0  =  50  m/s,  Acoustic  Interference  Model 


r/d 

Figure  56.  Radial  Density  Profile,  z/d=l,  U0  =  50  m/s,  Acoustic  Interference  Model 
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Figure  57.  Radial  Density  Profile,  z/d=2,  U0  =  50  m/s,  Acoustic  Interference  Model 


r/d 

Figure  58.  Radial  Density  Profile  Comparison,  UG  =  50  m/s,  1.0  ns,  Acoustic  Free  (AF) 

Model  and  Acoustic  Interference  Model 
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For  the  velocity  ratio  U„/Uj  =  6,  there  are  similar  results.  The  density  ratio  profiles  for 
the  acoustic  free  case  seen  in  Figures  59  to  64  and  the  acoustic  interference  models  in 
Figures  65  to  70  seem  to  follow  the  same  density  progression  seen  in  the  U0  =  0  m/s 
cases.  The  acoustic  free  density  ratios  for  U0  =  300  m/s  are  slightly  higher  than  those  for 
the  U0  =  0  m/s  cases.  The  density  ratios  for  the  acoustic  interference  model  also  favors 
the  U0  =  300  m/s  cases,  as  these  cases  see  higher  density  ratios  than  the  U0  =  0  m/s  cases. 
For  the  U0  =  0  m/s  case,  the  acoustic  interference  jet  is  wider  than  the  acoustic  free  jet 
and  the  data  receives  some  validation  from  Figure  51.  Figures  71,  72,  and  73  show  the 
acoustic  interference  jet  for  U0  =  300  m/s  is  not  as  wide  as  the  acoustic  free  jet.  The 
acoustic  free  jet  has  a  slightly  wider  radial  spread  than  does  the  acoustic  interference  jet. 
No  difference  is  apparent  in  the  axial  issuance  of  the  jet  into  the  chamber  for  either 
model.  The  radial  data  is  not  a  great  departure  from  previous  experimentation.  The 
acoustic  free  cases  present  a  slightly  wider  jet  in  the  radial  axis  with  little  difference  in 
the  axial  length  of  the  jet.  Figure  74  shows  shadowgraphs  from  the  experimental 
research  operating  at  higher  chamber  pressures  and  higher  velocity  ratios  and  small 
differences  can  be  seen  in  the  radial  width  of  the  acoustic  free  and  acoustic  interference 
jets  while  differences  in  the  axial  lengths  of  the  jets  is  fairly  indiscernible. 
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Figure  59.  Radial  Density  Profile  at  1.0  ns,  U0  =  300  m/s,  Acoustic  Free  Model 
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Figure  60.  Radial  Density  Profile  at  2.0  ns,  U0  =  300  m/s,  Acoustic  Free  Model 
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Figure  61.  Radial  Density  Profile  at  3.0  ns,  U0  =  300  m/s,  Acoustic  Free  Model 
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Figure  62.  Radial  Density  Profile  at  z/d=l,  Uo=300  m/s,  Acoustic  Free  Case 
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gure  63.  Radial  Density  Profile  at  z/d=2,  Uo=300  m/s,  Acoustic  Free  Case 


r/d 

Figure  64.  Radial  Density  Profile  at  z/d=4,  Uo=300  m/s,  Acoustic  Free  Case 
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Figure  66.  Radial  Density  Profile  at  2.0  ns,  U0  =  300  m/s,  Acoustic  Interference  Model 
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Figure  67.  Radial  Density  Profile  at  3.0  ns,  U0  =  300  m/s,  Acoustic  Interference  Model 


Figure  68.  Radial  Density  Profile  at  z/d=l,  Uo=300  m/s,  Acoustic  Interference  Model 
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69.  Radial  Density  Profile  at  z/d=2,  Uo=300  m/s,  Acoustic  Interference  Model 
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Figure  70.  Radial  Density  Profile  at  z/d=4,  Uo=300  m/s,  Acoustic  Interference  Model 
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Figure  71.  Radial  Density  Comparison,  U0  =  300  m/s,  1.0  ns,  Acoustic  Free  (AF)  Model 

and  Acoustic  Interference  Model 


r/d 

Figure  72.  Radial  Density  Comparison,  U0  =  300  m/s,  2.0  ns,  Acoustic  Free  (AF)  Model 

and  Acoustic  Interference  Model 
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Figure  74.  Shadowgraph  images  of  coaxial  jet  at  nearcritical  chamber  pressure 
(~3.5MPa)  and  at  the  high  outer-jet  temperature  (~  190K)  corresponding  to  cases  5-8 
(U0/Ui  =  ~  2.5,  ~7,  — 1 1,  and  -13.5).  The  acoustic  driver  is  off  for  images  in  the  top  row 
and  on  for  the  bottom  row  at  -3  KHz.  The  velocity  ratio  increases  from  left  to  right. 

(Davis,  2006:96) 
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V.  Conclusions  and  Recommendations 


Conclusions 

The  molecular  dynamic  simulation  tool  used  for  this  research  initially  studied  single  jet 
injector  flow.  This  research  required  changes  be  made  to  the  simulation  code.  Coaxial 
flow  and  acoustic  interference  were  incorporated  into  the  code.  The  new  code  has 
provided  simulation  data  for  single  and  coaxial  jet  flows,  with  and  without  acoustic 
interference  comparable  to  previous  experimental  data.  Multiple  velocity  ratios  were 
examined  and  their  results  were  compared  with  previous  experimental  data. 

The  original  coaxial  geometry  was  modified  to  provide  geometry  able  to  run  within  the 
present  computer  resource  limitations.  The  coaxial  geometry  provided  some  valuable 
data,  but  there  were  some  trends  seen  in  velocity  and  density  profiles  that  suggest  the 
geometry  may  not  be  ideal.  The  velocity  profiles  showed  acceleration  in  the  coaxial 
inner  and  outer  tube  as  the  flow  began.  The  coaxial  velocity  profiles  did  not  match  the 
theoretical  velocity  curves.  The  outer  velocity  curves  for  acoustic  free  and  acoustic 
interference  models  did  have  the  same  shape,  pointing  to  the  coaxial  tube  geometry  as  the 
reason  for  the  velocity  data  curves  obtained.  Acceleration  was  also  seen  in  the  single  jet 
case  when  the  acoustic  interference  was  applied.  The  effects  acoustic  interference  exert 
on  the  simulated  coaxial  jet  are  most  pronounced  in  the  mass  fraction  charts,  a  strong 
indication  of  mixing. 

Of  the  two  acoustic  interference  models  available  in  the  molecular  dynamic  tool,  only  the 
force  model  provided  expected  results.  The  discrete  particle  model  could  not  provide  any 
interference  due  to  the  size  of  the  simulations.  The  model  can  not  be  truly  validated 
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unless  simulation  size  is  increased  in  future  efforts.  Each  of  the  three  outer  to  inner  jet 
velocity  ratios  studied  provided  initial  data  used  to  validate  the  molecular  dynamic 
simulations.  The  single  jet  case  showed  less  mixing  when  acoustic  interference  was 
present.  The  acoustic  interference  case  for  the  single  jet  trailed  behind  the  mixing  of  the 
acoustic  free  case.  The  results  for  the  cases  with  a  coaxial  component  showed  improved 
mixing  when  the  acoustic  interference  was  applied.  Chamber  mixing  for  acoustic 
interference  models  exceeded  those  of  the  acoustic  free  case.  Compared  to  previous 
experimental  data,  there  was  correlation  showing  qualitative  agreement.  The  radial 
density  profiles  showed  results  that  matched  experimental  data.  Density  curves  for 
acoustic  interference  cases  showed  less  steep  declines  when  compared  to  acoustic  free 
cases.  The  research  is  in  the  early  stages  and  shows  progress  in  the  goal  of  providing 
simulations  modeling  desired  coaxial  and  acoustic  effects. 

Recommendations  for  Future  Research 

Future  endeavors  into  molecular  dynamic  simulation  with  the  existing  code  should  lead 
to  increased  testing  of  the  coaxial  jet  with  acoustic  effects.  More  variations  to  the  coaxial 
jet  ratio  of  D0/D;,  as  well  as  to  the  velocity  ratios  for  the  outer  and  inner  jets  need  to  be 
considered.  Presently,  the  simulation  takes  advantage  of  the  simplicity  of  modeling 
argon  and  nitrogen.  The  next  step  should  include  the  simulation  of  more  complex  fluids 
through  the  injectors.  As  the  simulations  become  more  involved,  the  required  computing 
capability  will  increase.  The  simulation  runs  for  this  research  pushed  the  computing 
envelope  and  required  modifications  in  order  to  be  completed  successfully.  The 
simulation  code  will  have  to  be  modified  to  create  an  environment  maximizing  parallel 
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computing  efficiency.  Simplified  calculation  algorithms  may  also  aid  in  decreasing 
calculation  times  for  simulations.  The  use  of  shared  memory  machines  will  reduce  the 
wall  run  times  by  50%  using  the  current  code. 

The  accelerations  out  of  the  tube  into  the  chamber  are  a  nano-hydrodynamic  effect 
coupled  with  the  mass  addition  subroutine  in  the  code.  There  is  a  net  pressure  caused  by 
the  size  of  the  tube.  These  problems  must  be  addressed  in  future  code  modifications. 
Possible  random  acoustic  interference  along  with  mechanics  to  deal  with  possible 
acoustic  reflections  from  chamber  walls  could  add  much  to  the  value  of  the  molecular 
dynamic  code. 


93 


Bibliography 


1.  Allen,  M.P.  and  D.J.  Tildesley.  Computer  Simulation  of  Liquid.  New  York: 
Oxford  University  Press.,  1987 

2.  Andersen,  Hans  C.  “Rattle:  A  “Velocity”  Version  of  the  Shake  Algorithm  for 
Molecular  Dynamics  Calculations,”  Journal  of  Computational  Physics,  52:  24-34 
(1983) 

3.  Au,  H.  and  N.W.M  Ko.  “Coaxial  Jets  of  Different  Mean  Velocity  Ratios,  Part  2,” 
Journal  of  Sound  and  Vibration,  1 16:  427-443  (1987) 

4.  Branam,  Richard  D.  Molecular  Dynamics  Simulation  of  Supercritical  Fluids. 

PhD.  Dissertation.  The  Pennsylvania  State  University,  University  Park  PA, 

August  2005 

5.  Branam,  R.,  G.  Schneider,  J.  Telaar,  and  W.  Mayer.  Length  Scales  in  Cryogenic 
Injection  at  Supercritical  Pressure.  Experiments  in  Fluids,  33:22-428  (2002) 

6.  Champagne,  F.  H.  and  I.  J.  Wygnanski.  “An  Experimental  Investigation  of 
Coaxial  Turbulent  Jet,”  Int.  Journal  of  Heat  Mass  Transfer,  14:  1445-1464 
(1971) 

7.  Danforth-Hanford,  Amanda,  P.  D.  O’Connor,  L.  N.  Long,  and  J.  B.  Anderson. 
Molecular  Relaxation  Simulations  in  Nonlinear  Acoustics  using  Direct  Simulation 
Monte  Carlo.  17th  International  Symposium  of  Nonlinear  Acoustics,  State 
College,  PA,  July  2005. 

8.  Danforth,  Amanda  and  Lyle  N.  Long.  “Nonlinear  acoustic  simulations  using 
direct  simulation  Monte  Carlo,”  Journal  of  Acoustic  Society  of  America,  1 16: 
1948-1955  (2004) 

9.  Davis,  Dustin  Wayne.  On  the  Behavior  of  a  Shear-Coaxial  Jet,  Spanning  Sub-  to 
Supercritical  Pressures,  with  and  without  an  Externally  Imposed  Transverse 
Acoustic  Field.  PhD.  Dissertation.  The  Pennsylvania  State  University,  University 
Park  PA,  May  2006 

10.  Davis,  Dustin  W.  and  Bruce  Chehroudi.  Measurements  in  an  Acoustically-Driven 
Coaxial  Jet  under  Sub-,  Near,  and  Supercritical  Conditions.  Engineering 
Research  Consultants,  Inc.,  Edwards  AFB  CA, 

1 1 .  Davis,  Dustin  and  Bruce  Chehroudi.  The  Effects  of  Pressure  and  Acoustic  Field 
on  a  Cryogenic  Coaxial  Jet.  Engineering  Research  Consultants,  Inc.,  Edwards 
AFB  CA, 


94 


12.  Haile,  J.M.  Molecular  Dynamics  Simulation:  Elementary  Methods.  John  Wiley 
&  Sons,  Inc.,  1992 

13.  Hadjiconstantinou,  Nicolas  G.  and  Alejandro  L.  Garcia.  “Molecular  simulations 
of  sound  wave  propagation  in  simple  gases,”  Physics  of  Fluids,  13:  1040-1046 
(2001) 

14.  Kamiadakis,  George,  Ali  Beskok,  Nayaran  Aluru,  Microflows  and  Nanoflows: 
Fundamentals  and  Simulation.,  Springer  Science+Business  Media,  Inc.,  2005 

15.  Ko,  N.  W.  M.  and  H  Au.  “Coaxial  Jets  of  Different  Mean  Velocity  Ratios,” 
Journal  of  Sound  and  Vibration,  100:  211-232  (1985) 

16.  Ludwig,  A.  “Molecular  Model  and  Basic  Physical  Variables”  Excerpt  from 
unpublished  article,  n.  pag. 

http://www.silcom.com/~aludwig/Physics/Main/Molecular  model.html 

17.  Ludwig,  A.  “Basic  Differential  Equations”  Excerpt  from  unpublished  article,  n. 
pag.  http://www.silcom.com/~aludwig/Physics/Main/Differential  equations.html 

18.  Ludwig,  A.  “Acoustic  Equations”  Excerpt  from  unpublished  article,  n.  pag. 
http://www.silcom.com/~aludwig/Physics/Main/Acoustic  equations.html 

19.  Ludwig,  A.  “Plane  Waves”  Excerpt  from  unpublished  article,  n.  pag. 
http://www.silcom.com/~aludwig/Physics/Main/Plane  waves.html 

20.  Ludwig,  A.  “Coupling  Power  out  of  a  Plane  Wave”  Excerpt  from  unpublished 
article,  n.  pag. 

http://www.silcom.com/~ahadwig/Physics/Main/Coupling  power.html 

21.  Mayer,  W.,  J.  Tellar,  R.  Branam,  and  G.  Schneider.  Characterization  of 
Cryogenic  Injection  at  Supercritical  Pressure.  AIAA  2001-3275,  37th 
AIAA/ASME/ASEE  Joint  Propulsion  Conference  and  Exhibit,  Salt  Lake  City, 
UT,  July  2001 

22.  Morton,  B.R.  “Coaxial  Turbulent  Jets,”  International  Journal  of  Heat  Mass 
Transfer,  5:  955-965  (1962) 

23.  Rapaport,  D.C.,  The  Art  of  Molecular  Dynamics  Simulation.  Cambridge 
University  Press,  1995 

24.  Schetz,  J.  A.,  Injection  and  Mixing  in  Turbulent  Flow.  New  York:  American 
Institute  of  Aeronautics  and  Astronautics,  1980. 


95 


25.  Strakey,  P.A.,  D.G.  Talley  and  J.J.  Hutt.  “Mixing  Characteristics  of  Coaxial 
Injectors  at  High  Gas/Liquid  Momentum  Ratios,”  Journal  of  Propulsion  and 
Power,  17:  402-410  (2001) 

26.  Villermaux,  Emmanuel.  “Mixing  and  Spray  Formation  in  Coaxial  Jets,”  Journal 
of  Propulsion  and  Power,  14:  807-816  (1998) 

27.  White,  Frank  M.,  Viscous  Fluid  Flow.,  McGraw-Hill,  Inc.,  1974 

28.  Wilcox,  David  C.,  Turbulent  Modeling  for  CFD.,  DCW  Industries.,  2006 

29.  Wochner,  Mark  S.,  Anthony  A.  Atchley  and  Victor  W.  Sparrow.  “Numerical 
simulation  of  finite  amplitude  wave  propagation  in  air  using  a  realistic 
atmospheric  absorption  model,”  Journal  of  Acoustic  Society  of  America,  118: 
2891-2898  (2005) 

30.  “Decibel  Scales:  Sound  Pressure  Fevel  and  Sound  Intensity  Fevel”  Excerpt  from 
unpublished  article,  n.  pag.  http://physics.mtsu.edu/~wmr/log  3  .htm 


96 


Vita 


Captain  Jermaine  Sailsman  graduated  from  Westtown  High  School  in  Westtown, 
Pennsylvania.  He  entered  undergraduate  studies  at  the  Georgia  Institute  of  Technology 
in  Atlanta,  Georgia  where  he  graduated  with  a  Bachelor  of  Science  degree  in  Aerospace 
Engineering  in  May  2000.  He  received  his  commissioned  through  Officer  Training 
School  in  Montgomery,  Alabama  in  August  2002. 

His  first  assignment  was  as  an  aircraft  data  engineer  at  the  National  Air  and  Space 
Intelligence  Center  at  Wright-Patterson  AFB.  In  August  2005,  he  completed  his  Master 
in  Business  Administration  from  Wright  State  University  with  a  dual  concentration  in 
Marketing  and  International  Business.  In  August  2005,  he  entered  the  Graduate  School 
of  Engineering  and  Management,  Air  Force  Institute  of  Technology.  Upon  graduation, 
he  will  be  assigned  to  Air  Force  Research  Fab,  Kirtland  AFB. 


97 


REPORT  DOCUMENTATION  PAGE 


Form  Approved 
OMB  No.  074-0188 


The  public  reporting  burden  for  this  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data 
sources,  gathering  and  maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any  other 
aspect  of  the  collection  of  information,  including  suggestions  for  reducing  this  burden  to  Department  of  Defense,  Washington  Headquarters  Services,  Directorate  for  Information 
Operations  and  Reports  (0704-0188),  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington,  VA  22202-4302.  Respondents  should  be  aware  that  notwithstanding  any  other 
provision  of  law,  no  person  shall  be  subject  to  an  penalty  for  failing  to  comply  with  a  collection  of  information  if  it  does  not  display  a  currently  valid  OMB  control  number. 
PLEASE  DO  NOT  RETURN  YOUR  FORM  TO  THE  ABOVE  ADDRESS. 


1 .  REPORT  DATE  (DD-MM-YYYY)  2.  REPORT  TYPE  3.  DATES  COVERED  (From  -  To) 

22-03-2007  Master’s  Thesis  August  2005  -  March  2007 


TITLE  AND  SUBTITLE 


Modeling  Acoustic  Effects  on  Shear  Coaxial  Jet  Flow 
Utilizing  Molecular  Dynamic  Simulation 


5a.  CONTRACT  NUMBER 

5b.  GRANT  NUMBER 

5c.  PROGRAM  ELEMENT  NUMBER 


6.  AUTHOR(S) 


Sailsman,  Jermaine  S.,  Captain,  USAF 


5d.  PROJECT  NUMBER 


5e.  TASK  NUMBER 


5f.  WORK  UNIT  NUMBER 


7.  PERFORMING  ORGANIZATION  NAMES(S)  AND  ADDRESS(S) 

Air  Force  Institute  of  Technology 

Graduate  School  of  Engineering  and  Management  (AFIT/EN) 
2950  Hobson  Way 
WPAFB  OH  45433-7765 


9.  SPONSORING/MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 

AFRL/AFOSR 

Attn:  Dr.  Mitat  Birkan 

875  North  Randolph  Street 

Arlington,  VA  22203  Commercial:703-696-7234 


12.  DISTRIBUTION/AVAILABILITY  STATEMENT 


8.  PERFORMING  ORGANIZATION 
REPORT  NUMBER 

AFIT/GA/EN  Y/07-M 1 6 


10.  SPONSOR/MONITOR’S 
ACRONYM(S) 


11.  SPONSOR/MONITOR’S  REPORT 
NUMBER(S) 


APPROVED  FOR  PUBLIC  RELEASE;  DISTRIBUTION  UNLIMITED. 


13.  SUPPLEMENTARY  NOTES 


14.  ABSTRACT 

The  purpose  of  this  research  is  to  determine  if  acoustical  effects  on  a  coaxial  shear  injection  jet  flow  can  be  modeled  through  the 
use  of  molecular  dynamic  simulation.  Molecular  dynamic  simulations  model  flows  as  a  group  of  interacting  particles.  The  flow 
in  this  research  was  simulated  using  nitrogen  molecules.  The  initial  task  involved  achieving  effective  geometry  for  a  simulated 
coaxial  jet.  The  coaxial  jet  geometry  was  driven  by  the  desire  for  simulations  to  operate  in  the  continuum  regime,  which  requires 
very  low  Knudsen  numbers.  Three  outer  to  inner  jet  ratios  of  0.0,  1.0,  and  6.0  were  examined  with  the  inner  jet  velocity 
maintained  constant  at  50  m/s. 

Velocity  profiles  in  the  coaxial  component  need  to  be  controlled  in  order  to  validate  the  continuum  flow.  Acoustic  interference  is 
introduced  into  the  simulation,  and  mixing  and  density  profiles  provide  valuable  information  into  the  how  the  flow  is  affected  by 
the  acoustic  interference.  Radial  density  profiles  also  provide  information  about  the  shape  the  jet  with  and  without  acoustic 
interference  as  it  exits  the  injector.  The  affects  of  acoustic  interference  for  most  cases  showed  good  agreement  with  the  previous 
experimental  data.  Results  showed  good  validation  of  the  simulation  and  warrants  more  in-depth  study. _ 


15.  SUBJECT  TERMS 

Acoustics,  Fluid  Flow,  Molecular  Dynamics,  Coaxial,  Interference 


98 


